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ABSTRACT 


Porous media are interconnected networks of void 
space, of multitude shapes and sizes. Attempts in developing 
an analytical expression for the fluid flow behaviour through 
them have not been satisfactory to-date. The classical model of 
the bundle of tubes is too simplified a model to be realistic 
and useful. In the present work the flow behaviour of porous 
media has been modelled by a network of interconnected tubep 
of different sizes. An attempt has been made to exfend the 
concept of relating pore-size distribution (both volumetric 
and number) to the capillary character of the network. The 
validity of ’ink-bottle’ effect has been tested by the 
comparison of capillary pressure data with the actual tube size 
distribution of the network. 

The extent of inter connection between the pores in a 
network is found to be a very important factor influencing the 
shape of the P^, curves. Therefore, any model of porous media 
based on bundle of tubes concept is likely to give erroneous 
results, Patt’s model, however, may be used with modified value 
of the exponent (between -2,5 to -3.5) to give better result 
for the number distribution. Similarly while correcting the 
measured pore-size distribution (volumetric) by Meyer’s method, 
it has been found that the relation Vo(].r^*^ improves the 
interpretation. These values for exponents have been obtained 
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for networks with- six tubes (on the average) meeting 
at a node point. The modified values of exponents perhaps 
correct for inter_connection between the pores. 




CHAPTER 1 


IHTROPUCTIOH \ 

, Porous media, whetlier in the form of natural rocks or 
particulate aggregates such as powder-packs, have quite a few 
interesting properties compared to the properties of their 
bulk solid phase. One such property, the flow of fluids through 
them, has been engaging the attention of scientists in many 
disciplines. The petroleum engineers are interested in knowing 
the quantity and the maximum possible rate of fluid withdrawal 
from reservoir rocks. The soil-scientists concern themselves 
with the distribution and movement of soil-moisture* Chemical 
engineers, likewise, dealing with the problem of drying, are 
faced with the question of determining the rate at which the 
moisture from the pores can move up to the surface. Similarly 
while analyzing the performance of a packed bed, it is pertinent 
to evaluate the effectiveness of mass-transfer as well as liquid 
hold-up. Internal void structure is important in deteimining 
the effectiveness of batteries as w'oll as the moisture movement 
capacity of the concrete and bricks. 

The total internal surface area, porosity or permeability 
do not provide adequate informa,tion about the complex pore 
structure of porous media. Quite often, it is necessary not 
only to know these cumulative properties but also the distri- 
bution of pore sizes. The shapes of the relative permeability 
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and. capillary pressure curves are riarkedly affected by pore- 
size distribution. Of the two, relative permeability does 
not seen to be suitable for detemining the pore-size distri- 
bution because of the averaging effect associated with the 
flow phenomena in the whole spectrum of pore-sizes. Capillary 
pressure, on the other hand, is a unique function of pore-size, 
interfacial tension and wettability. Keeping the latter two 
quantities constant it is possible to estimate the sizes and 
proportions of the contributing X-'ores from the successive higher 
values of capillary pressure. 

Capillary phenomena comes into plo.y v/henever a fluid/ 
fluid interface is in contact with a solid surface. Due to 
unequal interfacial tensions, fluid adjusts itself to attain 
theimodynanic equilibrium, thereby minimizing surface free 
energy. The well known Laplaco oquntion which relates the 
various quantities affecting capillary pressure is 



where = Capillary pressure, 

V' = Interfacial tension 
r^&r2= Two principal radii of curvature of the 
surface at any point. 

The implicit assumption that the sum of reciprocals 
of the two principal radii of curvature (mean curvature) is 
the same at every point on the equilibrium STirface, is a 
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mathematic ai consequence of the fact that the surface area, 
tends to a minimum (i.e* spherical or nearly so in most cases). 
All methods of pore-size determination, that make use of 
capillary pressure are based on the hypothesis that the mean 
curvature of an equilibrium liquid surface can somehow be 
related to the pore-size. This relationship is difficult to 
establish for any shapes other than that with the simplest 
geometry, for straight cylindrical capillaries, Laplace’s 


equation can be written as 

2 cos & 


( 1 - 2 ) 


c r 

where 0 is the angle of contact of the interface through the 
wetting phase and r is the radius of the capillary. This 
equation can correctly predict the poro-sizes and their distri- 
bution if bundle of parallel tubes is a valid model for porous 
media, for real system, this is an obvious over simplification 
and the equation will relate the capillary pressure to the 
radius of the largest opening leading to the pore rather than 
the radius of the pore itself. But in the absence of a better 
known method this pore entry radius is taken as a measure of 
the actual pore-size. Implied in this is the assumption that 
the actual pore-size is some unique function of the pore entry 
size for pores of all sizes. This however is difficult to 
justify. 

In capillary pressure experiment with porous media, 
the sample is gradually desaturated of the wetting phase by 
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injecting some non-wetting fluid (making known contact angle 


with the solid surface) at successive higher pressure. At 

equilibrium, the capillary pressure opposing the entry of 

non-wetting phase, P | is balanced by the externally applied 

pressure. This when put in Eq.(1.2) gives the radius of 

the smallest pore, r, inva,ded by the non-wetting phase. Knowing 

the saturation of the non-wetting phase in the sample at the 

stage, we arrive at the value of the fractional pore volume 

contained in the pores with ra.dius greater than or equal to r. 

By successively increasing P^, non-wetting phase is forced into 

smaller and smaller pores and the fractional pore volume 

corresponding to different values of r can be obtained. This 

method' has been universally used for the calculation of volumetric 

pore-size distribution in porous media. 

For fluid flow calculations, the pore number distribution 

is the desired quantity. To arrive at the number distribution 

from volume distribution, it is essential that ?;e make some 

assumption about the pore-geometry and its dim.ensions. For 

parallel tube model of porous media, the quantity in question 

is the relation between the tube radius and length. Based, on 

the assumption that 1 = c/-- , where C is proportionality constant 

1 

and r is radius, Fatt derives the number distribution function 

as follows? Sp 

J P ds 

o-j u 





(1.3) 
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where f(r 2 -r^) = frequenc.y of occurrence of pores in file 

interval r 2 -r^ 

S = saturation of the non-wetting phase 

There is considerable disagreement between different 

workers in the field about the value that ck. should assume, 

1 2 
¥/hereas Patt proposes a value of -1, Meyer implicitly takes 

this as 1 . Dallavalle^ seem,s to have some experimental results 

to back his suggestion of the value between 0 .and 1. In the 

absence of a conclusive evidence, this factor still remains 

uncerta,in. 

The method just discussed is essentially empirical in 
nature, for quantitative description of flov/ phenomena in 
porous media, a functional relation between the pore-size 
distribution and capillary pressure curve is needed. It will 
be our endeavour in the present work to exanine the existing 
models and suggest one which can reasonably interpret the 
experim.enta.1 vs saturation curves in terms of pore-size 
distribution. 


***** 



CHAPTER 2 


LITERATURE REVIEW MD AMALYSIS 

Hat-ural porous systems can be viewed as inter-connected 
network of voids of different shapes and sizes separated by 
interconnecting links or necks. Often the behaviour of such 
a system is approximated by a bundle of capillaries^. Duliien 

5 

and Batra have collected a comprehensive list of references 
of the attempts made in this direction in their review article. 

It is apparent that very fev; real porous media will contain 
pores of the straight cylindrical shapes, and therefore the 
pore radii calculated from the capillary pressure have to be 
considered as equivalent pore radii. This is defined as the 
radius of a straight cylindrical capillary that would give rise 
to the same capillary pressure as the measured value. Since a 
variety of different capillary shapes may give rise to the same 
capillary pressure, the pore dimensions calculated from capillary 
pressure measurements are at best semi-quantitative. The equation 
relating properties of porous media to the radius distribution 
of the equivalent bundle of tubes have been given by a number 
of authors . However, apart from the mathematical simplicity 
associated, the most glaring weakness of the model is the 
prediction of highly anisotropic properties. Porous systems 
on the other hand are known to be isotropic in nature. Inspite 
of its success in correlating certain properties viz, permeability 
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formation factor etc., the model has failed to approximate 
the capillary pressure vs saturation and its hysteresis, 

For more complex systems, such as random packing of 
non-spherical particles, Carman*^ suggests that tv/ice the 
hydraulic radius is the suitable value of r in the following 
expression for P ; 


2 Y cos Q 
r 


( 2 . 1 ) 


where the hydraulic radius is defined as the ratio of porosity 
f to the surface area A in unit bulk volume. A more accurate 
account of the curvature for different geometry of the capillary 
is given by Gregg and Sing . On substitution Bq..(2,l) is 
reduced to . . 


^ .. . . yt c os . S _A ( 2 . 2 ) 

Substituting the expression for the specific surface area 

1 1 

in terms of permeability K and porosity, from the Kozeny 
equation, Eq.(2,2) becomes 


_ (J.)V 


1 i 


(2.3) 


where the threshold pressure P^ has been substituted for P^. 
k is an empirical textural constant of the medium and accounts 
for the deviation of a porous system from the model. 

This was obviously a gross simplification, and Carman 
himself points out that, for porous bodies with widely varying 
pore-size distribution, calculations based on Kozeny-Carman 
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model can be quite wrong. 


Sphere Pack Model s 


The earliest attempt in approximating the flov; behaviour 

of porous media was their simulation by packs of uniform spheres. 

But the complexity of the pore-geometry in sphere-packs prevented 

the derivation of an accurate and meaningful description of 

flow behaviour through them. Recognizing the difficulty, 

1 1 

Koseny inductively developed an equation relating permeability 
to porosity and internal surface area of a sphere-pack. 


where 


K = 


C 

T? 


f = porosity 


(2.4) 


T = tortuosity 
S = specific surface area 

12 

Kozeny equation, as modified by Carman and given as 


K 


k S (l-f)‘^ 


(2.5) 


did not prove to be of much value, aside from its use in 

estimation of the surface area of powders. 

A comprehensive review of the use of sphere -packs to 

model the flow behaviour of packed beds, has been given by 

1 5 

Haughey and Beveridge , In their analysis, the emphasis has 
been in arriving at the mean fractional voidage for both 
regular and random packing, on the basis of the position coordi- 
nates and the coordination numbers, rather than the pore « 
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geometry. Generally speaking, particle shape and. size 
distribution are the two factors most likely to affect the 
packing structure and its properties. Some work has been 
carried out on the densest packing of non-spherical particles. 
But the multiplicity of packing types with the same coordination 
number and vice versa, a,ccompanied by the size distribution of 
particles further complicate the pore structure. Hence the 
task of analytical description is rendered almost hopeless. 

The problem of capillary hysteresis in uniform sphere 

i 4* 15 

packs was studied by Kruyer Prevel and Kressley Meyer and 

1 6 1 T 18 

Stowe ’ and more recently, by Melrose . Their analysis 

points out that the curvature of the invading interface of 

non-wetting phase is generally larger resulting in higher value 

of capillary pressure compared to the retracting interface 

associated with imbibition.' This, according to them, may be 

the reason for the porous media exhibiting capillary hysteresis, 

IQ 

Haar and Wygal " view any porous media as a random 
mixture of spherical particles and propose a model based on the 
principle of averaging the properties of the basic sets to 
arrive at the mixture-properties. The problem essentially is 
finding of appropriate weighting factors to average the 
basic set properties, determined experimentally. Expressions 
for porosity, permeability and capillary pressure behaviour 
of the porous media have been developed. But the experimental 
support to this model is too meagre to be accepted. 
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Network Model s 

The description of single and multiphase fluid flow 

through capillary pores in porous media has been unsatisfactory 

because the pore-geometry is too complex to be described 

20 

adequately by analytical expression. Markin used a model 
consisting of a system of randomly arranged intersecting voids 
with circular cross-section and continuously varying radius. 

No more than three branches could converge at any point in the 
netjA?ork, This model, however, does not seem to have gained 
the favour of the research workers in the field, 

Patt replaced each pore with a cylinder in his two 
dimensional network with different grid patterns. He could 
incorpois-te random assignment of parameters, such as pore-size 
distribution and number of cylinders meeting at a node point, 
in the model. Fluid flow in the network could be followed as 
one pore after another was desaturated. 

All of Fatt’s calculations were performed assuming 
a continuous wetting phase in the tubes. As a consequence, 
no wetting phase could be trapped in a de saturation process. 
Dodd and Kiel extended this v;ork by applying desaturation 
steps such that wetting fluid could be trapped. Cylindrical 
pores were assumed to contain only one fluid at a time; thus 
allowing the displaced fluid to be trapped if no continuous 
path to the effluent end was available. Thus at any time 
the pore could be either full of wetting phase or non-wetting 
phase. 
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Singhal modified fatt's model ty providing for tlie 
possilDilities of different flov^ reg3.mes in different channels 
during desaturation process. The four basic flow regimes viz, 
single phase flow, displacement, annular and slug flow were 
allowed to occur any where in the model whenever suitable 
conditions existed. In addition, he substituted triangular 
capillaries for circular shaped tubes to account for funicular 
flow regime over extended range of saturation. A trapping 
factor, evaluated from the comparison of the model with the 
experimental results, was used to effect trapping of fluid. 
This model has duplicated the relative permeability behaviour 
of porous media reasonably well and we propose to use the same 
model for approximating its capillary pressure vs, saturation 
behaviour. 

Other Model : 

2d 

Haring and Greenkom have proposed a parametric 
statistical model, from the theoretical considerations, to 
calculate the macroscopic properties such as saturation, 
permeability and dispersion coefficient in porous media. 

Porous medium is approximated by a large number of randomly 
oriented straight, cylindrical pores with randomly varying 
radius and length. The relationship between saturation and 
P_ contains two constants which are the parameters of the 
pore-size distribution function. They maintain that the two 
parameter incomplete Beta function can adequately represent 
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most of the pore-siae distributions encountered in porous 
media* Non-uniformity of a medium can be incorporated by 
adjusting the magnitude of the two parameters in the distri- 
bution function. Of all the models proposed to-date this 
appears to be the most promising for analytical description 
of flow phenomena as well as mass transfer through porous 
media, 

24 

Pay takes et al, while developing a model for the packed 
bed, point out some in-consistency in Haring and Greenkorn’s model. 
Whereas in Haring’s model the minimum tube size is implicit^ly 
assumed as zero, the smallest pore-size in a packed bed is 
definitely not zero. As a' consequence, the average pore diame^i^ 
calculated from this model are smaller than the smallest pore 
calculated from sphere pack model. The calculated pore length 
is about equal to the half the average grain diameter of the 
pack. Hov;cver, this minor inconsistency in the model can be 
removed by substituting for the finite value of the smallest 
pore-size for packed bed and thereby slightly modifying the 
model. 

In their own model, Paytakes et al assume convergent- 
divergent circular tubes as the basic units of which a packed 
bed is assumed to comprise of. These basic units do not have 
any lateral connection thus allowing fluid flow only in the 
vertical direction. By introducing convergent-divergent flow 
channels, they seem to have accounted for inertia effects in 
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flow path. However, this could not he a realistic approach 
for porous media where interconnections are quite important 
in considering the flow hehaviour. 

2.2 Pore-Size Distrihution ; 

It is obvious that pore-size distribution is very 

important in describing fluid-flow through porous media, 

25 

Ritter and Drake proposed a method of determining the pore- 
size distribution from mercury porosimetry data. The basis 
of this method is the concept that mercury is forced into 
gradually smaller pores against capillary forces as the pressure 
is increased. The fact, that such information does not directly 
give the pore-size distribution is obvious. However, pore-size 
distribution obtained from Nitrogen adsorption agrees well with 
that calculated from mercury porosimeter data on the assumption 
of parallel tube model. ° Mercury porosimetry will not detect 
the dimensions of those pores that are accessible only through 
necks narrower than the pore itself. The volume, of these pores 
of limited accessibility will be erroneously assigned the 
dimensions of the largest neck leading to these pores. They 
analyzed a number of natural substances a,nd concluded that the 
pore size distribution can be represented by three parameter 
modified Maxwellian distribution function. Haring’s two para- 
meter distribution function is simpler and should be preferred. 
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Meyer attempted to correct this distortion in the 
pore-size distribution due to the presence of the pores of 
Limited accessibilitj?'. On the assumption that the pore sizes 
follow Poisson distribution, he calculated the probability of 
a pore of radius r^ not being connected to the pores of radius 
r r^ . This probability was used to correct the mercury 
porosimeter values for the 'Ink-bottle' effect. The details 
of the derivation and the computational method is given in 
Appendix B, 

27 

DuHi.en proposed to determine 'two dimensional’ pore- 
size distribution by injecting ¥/ood's metal into the sample 
and using micrographic analysis. Apart from the experimental 
scheme, he has not published his results yet. 


2 , 3 Characterization of Porous Media on the Basis of Capillary 
Pressure Curve ; 

The earliest attempt in this direction was made by 

p Q p Q , 

Leverett ’ ^ who put forward a semi-empirical relation 




(I)* 


( 2 . 6 ) 


'w' ''Y cos 0 f 

v/here ” capillary pressure function and is dimension- 

less function of wretting phase saturation 
on the basis of dimensional analysis, Leverett 's data showed 
that a plot of J(S^) vs (wetting phase saturation) yielded 
a unique curve which described the capillary retention of wetting 
liquid in clean, unconsolidated sand. 


15 


30 

Rose and Bruce extended Kozeny-Carman treatment to 
relate capillary pressure with porosity, permeability and 
textural constant known as Kozeny constant. The expression 
suggested iss 


'Yoos 9 


(e) 


( 2 . 7 ) 


where P, 


T 


threshold pressiire i.e. lim P 

Sv4*tyl 


•T 


But these expressions, apart from classifying the porous 

media and in some cases predicting flow behaviour, did not 

give any insight in the distribution of the pore-sizes, 

P'5 

Recently Haring and Grcenkorn" have derived an expression for 
wetting phase saturation from P^ curve as 


1 

R * 

rc ^ .X+2 

> V 'P* ^ 


S. 


0^ 


.p*- 

c 


(1 - PV)'® 

. 0 c 


w 


[ ( 1 - 


(2.8) 


0 c 


where P* 




p 

^T 

threshold capillary pressure 
parameters of pore size distribution. 


Thus in case we ea?e able to estimate /3 from P^ vs 

curve, we really have an expression for pore-size distribution 


5 '(x), 
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I dx 

F(x*) = -ij — (2.9) 

J x”"' (l-x)^^ dx 
0 

where x* = ^/^max ’ ratio of pore-size to the largest 

pore size. 

This, F(x) is then used to calculate the wetting fluid 
saturation, for any by the Eq, (2.8). 

2.4 Plan of Attack s 

Prctn the analysis of tlie past work, it is obvious that 

the relation between the po le -structure, its distribution and 

the capillary characteristics of porous media is not ivell under- 
2 

stood, Meyer assumed Poisson distribution and suggested a 

method for applying correction fbr ’ink-bottle’ effect. Haring 
25 

and Greenkorn suggested a statistical model for capillary 

pressure which could be useful in mass transfer calculations. 

But they approximated the pores by straight cylindrical tubes 

beside assuming unimodal distribution of the poies. In real 

porous media these are not strictly valid assumptions. In the 

interpretation of P_ curve for pore-size distribution, Patt 

used an arbitrary factor i^,) which could take the value -1, 

There seems to be considerable disagreement between different 

2 

workers regarding the value of this factor, , Meyer had 
implicitly assumed the value 1 for in his calculation of 
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of probability function for the pores of limited accessibility. 

To find a method which can interpret P curves satis- 
factorily, we will do the follov;ing; 

a) networks will be constructed for different pore-size 
distributions and by allowing them to imbibe nonwetting fluid 
under successively higher pressure, P^ curves will be generated 
for different distributions, 

b) An attempt will be made to v/orkout a procedure for 

interpreting P in terms of the pore-size distribution. To 

arrive at the number distribution of the ' po re-sizes from P„ 

1 

curve, Patt ’ s model with different values of y/,. will be tried 
and the result willhe compared with the actual input to the 
network. 

c) Existence of ’ink-bottle' effect in a real porous 

2 

system v/ill be studied by the method proposed by Meyer . The 
corrected pore-size distribution (volumetric) will be matched 
with the one estimated from the P curve to determine whether 
the effect is significant in natural system. Corrected pore- 
sizo distribution from the network P^ curve wall be compared with 
the actual distribution to test the reliability of the approach. 
To this end, different expressions for characteristic volume, 
in Meyer's method, in tertifis of power to pore radius, will be 
tried and compared with the actual input to the network. 

d) To obtain the pore-size distribution from an actual 

23 

P„ curve. Haring and Greenkorn model will be used to estimate 
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tlie parameters of distribution by non-linear least square 
method. Network model will be used to determine how close 
the estimated distribution is to the actual distribution. 


***** 



CEigTER 3 


MORELS ESED IE THE ERESEET WORK 
3.1 Generation of Curve- from Known Pore-Size Distribution-; 

' """ '■’"“■'"■f'l'i -I '■'■'I - TT---I.ri,.i 0^ " IIUJI r--iri "i ti|-| - i '-r -'rriijr Tr--ii inr"rir -r- -iii i -fi'ii ;■ i tin i 'T- it- ' T-ir-ni, ir - tV'IT' 1 

Use of network model for generation of P. curve consists 
of two steps; 

a) Construction of a netvrork simulating the porous 

material with a knovm pore-size distribution. This idea was 

originally suggested by Patt and the methodology followed here 

22 

was developed by Singhal . In the network , pores could merge 
at and branchout from a node point, A network of the size 
6x8 was used in our simulation in which the nodes were arranged 
in square grid pattern. Thus at each node point a mfiximum'- of 
eight tubes cou-ld meet, Patt suggested, from experimental 

observations that each tube was connected to 10 other tubes 

onIKc. 

implying that 6 tubes could meet at a point. This fact 

was incorporated in the model by modifying the cumulative 
pore-size distribution obtained from capillary pressure curves, 
so that 25^ of the tubes have zero radius. This was achieved 
by compressing the real distribution between 0.25 to 1.0 instead 
of 0.0 to 1,0 by the following expression; 

1 ’ = 0.751 + 0.25 

1’ = new value of the distribution 
1 = original value 


where 
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Whereas the tube radius and location were assigned in 

manner by Monte-Carlo technique, the tube length could 

assume only two values viz. equal to the side length of the 

square or its diagonal length. The grid length was arbitrarily 

chosen. The list of the input to a typical network simulation 

program is given in the Table 1, The details of simula.tion 

22 

has been discussed by Singhal and for the sake of completion, 

the relevant points are included in the Appendix A, 

h) Generation of Ca:;'ill’ary Pressure Curves Capillary 

pressure characteristics of the network is evaluated by modifying 

22 

and adapting the computer program written by Singhal Some 
additional checks have been introduced in the program to verify 
that material balance for each fluid as well as for the total 
fluid is satisfied at each stage of computational cycle. 

Initially the network is assunied to be full of v/etting 
fluid and is allowed to be displaced. by non-wetting fluid, A 
fixed pressure drop ( JiP) is applied across the network and 
it is allowed to desaturste the wetting fluid in a quasi- 
static manner. This process continues till outflow from the 
network in a particular time step is negligibly small; implying 
equilibrium- fluid distribution in the model. Theoretically 
in a quasi-static process the time stop should be as small as 
possible. But this would mean larger computation time. Hence 
we used the time steps of the order of 200 /'.sec. initially 
and later increased it as the desaturation process proceeded. 
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OCCUY 

Instead of v>?aiting indefinitely for zero outfloY//,in a time 
interval, we cut short the calculation for a P whenever 
the outflo?; in a time step is less than or equal to 1?^ of the 
outflow in the first time step. The saturation of the wetting 
phase in the model is then calculated thus obtaining one point 
on the vs S„, curve. Pressure drop ( Zi P) across the model 
is then successively increased and in each case the equilibrium 
wetting phase saturation is obtained. In this way the whole 

P vs S curve can be obtained over any range of pressure and 

c* w 

saturation. 

Singhal^"" used a trapping factor (0,2/^i^^^) , where 
Xa. and are the viscosities of non-wetting and wetting 
fluid respectively, to trap the fluids in the pores during 
desaturation. This was merely a correction factor to match 
the model behaviour with actual samples. In static process 
like capillary desaturation, the viscosity of the fluid is 
not of any significance as far as equilibrium saturation is 
concerned. Equilibrium saturation of the wetting fluid depends 
solely on the curvature of the remaining fluid interfaces in 
various tubes. In our simulation of the network the tubes have 
assumed the size equivalent to the volumetric average radius 
of the pores. In real porous systems ( granular) the pores 
are convergent-divergent in shape. As a result the neck radius 
in any tube is definitely smaller than the volumetric average 
radius. As the capillary pressure of any tube is determined 
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by the size of the narrowest constriction in the tube, the ' 
actual capillary pressure will always be larger than that 
calculated from the average radius. 

In the present work, the trapping factor is the ratio 
of the volumetric average pore radius to , the neck radius for 
each tube, thus implicitly simulating the convergent - divergent 
shape of the pore channels. At every stage the pore-entry 
capillary pressure is modified by this factor to calculate 
at the neck. This pressure is then used to determine whether 
the interface will move through the pore or will be blocked, 

3,2 Correction for "ink-Bottle” Effect ; 

In natural porous systems pores are surrounded by 

2 

smaller pores; recognized as "ink-bottle*' effect by Meyer , 

This results in assigning larger volume to the smaller pores 
while the volumetric share of larger pores is less than the 
actual volume. Evidently the use of bundle of tubes model to 
calculate pore size distribution is likely to be erroneous, 
Meyer‘S attempted to correct the measured data for this effect 
by calculating the probability of a larger pore not being 
connected through equal or larger pores to the mercury source 
in mercury porosimetry. He assumed that pore-siees are 
randomly distributed and follow Poisson's distribution. In 
deriving the probability overlooks the fact 

that the pores are mutually exclusive and asserts that this 
will not introduce any serious erior in calculation since only 
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the 'zero’ term in Poisson's distribution is used for calculation. 
Also, assumed randomness of pore sizes exclude, from consideration 
clusters of very large or very small pores such as fractured 
or vuggy systom. 

Other simplifying assumptions made by Meyer are; 

(i) that the centre of an r i.X)re lying in volume v(r^,r) will 
ensure its being connected to an r^ pore, (ii) that only those 
pores which are isolated from all equal or larger pores will 
fail to fill and (iii) volume of a pore may be expressed as 
Kr"' where K is a constant. 

The first assumption is valid for spherical pores and 
seems intuitively sound for others. The implication of the 
second assumption is that isolated groups of mutually connected 
pores are not being considered. This, of course, is a serious 
simplification but necessary for handling the difficult subject 
of aggregates in a random distribution. It is hoped that the 
correction applied for individual isolated pores is the major 
correction needed. The third assumptions implies that larger 
pores have larger length a£5socio.ted with them. This sounds 
reasonable as the larger pores are associated with larger grains 
for uniform packing. This assumption has one snag that K is 
not a constant. However, K as a constant appears only under 
the integral sign and cancels out. The assumption, in fact, 
is not that K is, a constant but that its mean value in one 
integral is equal to that in the other. 
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Hence, lie pjcoposed that measured pore-size distrihution 
should be corrected by this method to obtain actual distri- 
bution for granular porous system. Details of the method is 
given in Appendix B. 

Computer program, written for this scheme, is given in 
Appendix D. 

3.3 Haring and Green horn Model ; 

23 

Haring and Grecnkom"' have suggested the use of two 
parameter incom.plete Beta function as the universal distribution 
function for the pore-size in -a porous body. It was claimed 
by them, and rightly so, that almost any unimodal distribution, 
whether symmetric or skew, could be approximated very closely 
by the adjustment of the two parameters. Normal distribution 
of the pore-size was not favoured for the following two reasons; 
i) the distribution is symmetrical, 
ii) the independent variable varies from -i;K5to 5 

none of which is strictly true for a porous media. 

One more implicit assumption in their approximation 
is that the ratio of minimum to maximum pore size is zero, 
which is not always true. That seems to he the reason, why 
for a granular pack the average pore size computed on the 
basis of this distribution function comes to even less than 
the minimum pore-size computed from geometrical considerations 
of granular packing. Also, the average pore length is calculated 
to be equal to half the average particle diameter; vfnile ideally 
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it should he equal. This fact has been pointed out by Paytakes, 
.24 

Tien and Turian . However, this inconsistency may be over 
come by slightly modifying the distribution function to ; 

-7*^ ft 

i (1-x)'^' dx 

g(r*) = ^ (5.1) 

l(x-x (l-x)^'^ dx 
% Tn 

where r* = dimensionless radius of the pore. The ratio 

of pore radius to the largest pore radius, 

= The ratio of smallest to the largest pore 
radii 

^ , ft = Parameters of distribution. 

Based on the assumption that the actual pore-size 

distribution can be approximated reasonably by incomplete Beta 

23 

function, Haring et al arrived at the following expressions 
for the vifetting fluid saturation, 


where 


1_ 

P* + 3)1 

^ " f ^ '( ^..+2)1 13 I 
5 


(1-x)^’^ dx 


= B 


( 1 . 

^ p* 

c 


.oC+2 A 



( 3 . 2 ) 


( 3 . 3 ) 


1 

B( w ? cic+'S j /i ) is incomplete Beta function with 
. i^c 

■y''^+2 and as parameters,. 

P* Dimensionless capillary pressure, 

Pj Threshold capillary pressure i.e. P^, when 


•™?' 1 , 
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In this work we shall try to fit the Eq.( 3 . 3 ) hy 
non-linear least square method. Standard- I BI program SDA- 
3094 will be used for this purpose. 





CHAPTER 4 


RESTJETS MJ) DISCHSSIOH 

As already mentioned in earlier cliapters network 
seems to He a good model for predicting the flow behaviour 
of porous r^iedia, A part from its closer and better physical 
similarity w'ith porous media, capillary pressure curves 
obtained from it look quite realistic. However, there may 
be some weaknesses in this model which will be discussed 
later in this chapter. 

Capillary pressure (P^^) vs wetting fluid saturation (S^) 

from the networks are given in Pigs. 2a to 2f. Each network 

was generated from a specified pore-size distribution, vdiich 

was assumed to represent a hypothetical porous medium. Pore- 

size distribution used to generate networks were; arbitrary 

distribution, Beta distribution and Poisson distribution. 

Besides, a couple of biraodal arbitrary distributions were also 

used. The netvrorks with bimodal distribution were used to find 

out if the P. curve could give some indication of the multi- 

modality in the distribution. In each case, except the first 

network, the largest pore-size used Vi?as 20 JX , The input 

distribution was compressed in each case between 0.25 to 1.0 

22 

to accommodate Patt's Beta factor as suggested by Singhal . 

Inspite of the very distinct and extreme forms of the 
distributions used to generate the networks, th« caioulated 
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ciir-ves were not significantly different from eacli other. 

Thus giwing the impression that capillary pressure is not a 
very strong function of the pore-size distribution in the net- 
work, Also threshold pressures (Pq,) obtained from the network 
P^ curves (corresponding to = 1.0) do not always agree with 
that predicted from the bundle of tubes model. Threshold 
pressure, from the definition, depends on the size of the 
largest tube in the input face of the netv/ork. But, in view 
of the fact that tubes are assigned randomly at different node 
points, it is possible that the largest size tube may not 
occur at the input face,, In that case Prp predicted from the 
bundle of tubes model will bo different from that evaluated 
from the network. Prom the inspection of Pig. I.e, it is 
obvious that the largest tube size (20,4A- ) does not appear in 
the input face and henco the calculated P^, from the network 
is higher than that predicted for 20 , 01 - tube radius (0.259 psi). 
Also, during the computation of P curve from networks, more 
than one fluid/fluid interfaces were observed in a continuous 
flow channel. This factor has not been accounted for in the 
bundle of tubes model and can be one of the reasons resulting 
in the difference in the two values of 

Other features of the network likely to have influenced 

the P curves are? (i) number of tubes (ii) ratio of the tube 
c 

length to the radius (iii) the extent of inter-connection 
between the tubesand (iv) existence of different flow regimes 
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in flow channels. 
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(i) fatt estimated number of pores in a 1/2” cube sample 
of sandstone to be of the order of 10,000. In our networks, 
the number of tubes varied between 89 to 101. ITiether such 
a small number of tubes (pores) could reasonably approximate 
the flow behaviour of real porous media with huge number of 
pores, is open to doubt, V/ith such a small number of tubes, 
it is the relative position of the tubes which matters 

more them the size distribution. If the investigation could 
bo repeated with larger sizes of networks (500 to 1000 tubes), 
the result may throw some light over thijii question. Our 
limitation in computer memory and time prevented us from 
trying bigger networks, Patt a.sserts that tubes more than 
200 in number do not substantially alfect the nature of 
curves from networks. Kovvever, he aid not compare any 
experimental P_ curve w'ith that of a network to prove his 
point that tube number is not a critical factor. 

(ii) The second factor, which seems significant, is the 
ratio of the pore length to the pore radius. In our networks, 
this ratio varies from o,bout 6 to 60 depending on the size 
and the location of the tubes relative to the grid. Earlier 
investigators suggest thef the length of a poie varies, accord- 
ing to some power of r,. This factor varie.s between 0 and 1, 
according to Dallavalle^, Patt”* puts this value at -1 , In 
one of our calculations grid size was reduced to about l/5th 
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of the original value thus reducing the maximum value’ of 
length to the radius ratio from 60 to 12, The two P. curves 
(Pig.lh) are not significantly different, Further investi- 
gation with different values of tl^e grid-size is desirable 
before- drawing any conclusion. 

(iii) The third factor to be considered is the Beta factor 

1 99 

suggested by Patt and used by Singhal in his simulation 
of networks. This factor essentially represents the extent 
of interconnection between the pores in porous media. Applied 
to networks, this implies average number of tubes meeting at 
a node point, Patt , from his experimental observations, points 
out that in a sandstone rock each pore is connected to 10 other 
pores thereby implying tjiat on the average 6 tubes should meet 
at a node point. In Pig* 2h the three P curves have been 
obtained with different values of the Beta factor. It is 
obvious that the two P curves, for Beta factors 6 5: 8, are 
not very different from each other, the P curve for Beta factor 
equal to 4 is significantly different from the other two. Hence 
it may be inferred that the extent of interconnection in a 
network, represented by Beta factor, is quite important and 
is likely to be different for different media. In our simulation 
this factor has been taken equal to 6. 

( iv) There is experimental evidence to support the existence 
of different flow regimes in flow paths of physical models of 
porous media where one fluid displaces the other. In such 
systems, the interplay of viscous and capillary forces comes 
into existence and the flow regimes depend on A ? and flow 
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rate. This phenomenon is accounted for in the network but 
not in the bundle of tubes model. This could be another weak- 
ness of the bundle of tubes model-. However, for mercury- 
vacuum system this situation, of different flow regimes in 
flow paths, is not likely to occur. 

Simulated networks used for generation of curves were 

2 

testeti for«2oodw2ss of -fldby test. The pore-size distributions 

of the simulated networks were not found to be similar to 

the input distributions within a reasonable confidence limit, 

V/e suspect that the scheme of pseudo random numbers used in the 

simulation might not have been truly random. Consequently it 

was decided to use the actual distributions of the networks 

rather than the input distributions v/hich they were supposed 

to represent, Any further v'ork in this direction should ensure 

2 

that the simulated network does pass "v' test to represent 
the input distribution. Hov/c-ver, the weaknesses mentioned above 
may not be very serious in view of the results obtained from 
the network model. 

Since our netw'orks have discrete pore-sizes, ideally 
some step type P^ curves would have been expected. But from 
the limited number of the points on the P^ vs plot, a smooth 
curve could be dravm. It appears that the pore-sizes may not 
be the exclusive factor determining the value of capillary 
pressure in a porous medium. Other factors like connectivity 
etc might also bo important, . 
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The idea of relating to the size of the pores has 
its root in the bundle of tubes model of porous media. Even 
the model proposed by Haring and Greenkorn^^ are the modifi- 
cation of the same model, providing random' orientation of the tu 
tubes and consequent inter-connections. Therefore to correct 
for inter-connection between the pores, some modification is 
called for in the interpretation of P. curves based on this 
model. 

The effect of "ink-bottle" phenomena seems to be signi- 
■ficant in the networks, resulting in the cumulative pore-size 
distribution (volumetric) obtained from P curves, being 
quite different form the actual distributions. Actual distri- 
butions were computed by counting the number of tubes of ^ 
different size and calculating their volumes. The difference 
in the two distributions can be clearly seen in the Pig, 3a to 
3f . Attempts were- made to correct the distributions (obtained 
from P^ curves) by Meyer's method for "ink-bottle" effect. 

In using this method, slight modification was made in the scheme 
given in the Appendix B, As is obvious, in a network, tube 
lengths are independent of radius. As a result the pore volume 
as well as the characteristic volume defined in Meyer’s method 
becomes proportional to the radius squared rather than radius 
cubed. The corrected distributions (volumetric) is compared 
with the actual and that obtained from P^ curves in Pigs, 3a 
to 3f , 
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Effect of the pores of limited accessihility ("ink- 
bottle" effect) on the volumetric distribution appears quite- 
pronoimced in some of the natural substances investigated by 
us. The substances included porous plate, fritted glass, pelleted 
gel etc,. Meyer's method as given in the Appendix B was used 
for these calculations and the result is given in Fig. 4a to 4di' 

25 

Pjj data for these substances were obtained from Ritter and Drake 

and are reproduced in the Table 3. The same method was applied to 

data of some sandstone samples obtained from literature but 

no appreciable difference in the experimental and corrected 

distributions was observed? implying that this phenomenon is not 

very significant in those samples. In an attempt to closely 

approximate the actual pore-size distribution, the exponent of 

r in Meyer's method was given a number of values (3 to 1,5)* The 

curves with -the exponent equal to 1,5 were observed to be better 

than others in approximating the actual curves (Figs. 5a to 5^')* 

0 ^ 

Physically this means that the pore length iscjCr - . This 
distortion from the expected relation of loCr° may be due to a 
number of factors, like some bias in our networks, and to a 
greater degree on the extent of inter connection between the pores. 
Attempts were also made to estimate number pore-size 
distribution from P„ vs curves obtained from the networks. 

In the absence of an analytical expression incorporating net- 
work concept of porous media, investigations were carried out 
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to find out if Patt’s model, of bundle of tubes, could be 
useful. In that model (Eq, 1.3) Patt"^ used an exponent which 

related the pore length to the radius. The same equation, with 
the different values of c>C ? was used to calculate pore-size 
distributions (number) from P curves of the networks. It is 

* L/ 

apparent from Pigs. 5a to 5f that for the values of lying 

between -2.5 to -3.5 (mostly -3.0) the calculated distributions 

are reasonably close to the actual distribution. The physical 

— 3 0 

implication of this value- is , that the pore length (l) r 
In our networks the tube lengths were independent of the radius. 
But because of the fact that we have distinctly two groups of 
lengths, side lengths and diagonals of the grid, in our net- 
works, pore lengths could not be said entirely independent of 
the radius. Also, the small sample size, extent of inter- 
connections and the non -randomness of the pseudo random numbers 
used in simulation, might have introduced some fictitious 
relation betvrcen 1 and r. 

23 

Haring and Greonkorn model relating pore-size distri- 
bution (number) to capillary pressure and saturation, was 
tested by estimating the parameters of the distribution from 
P^ vs S , curves by non-linear least square technique. These 
parameters ?;ere put back in incomplete Beta distribution 
function to calculate the number distribution. These calculated 
pore-size distributions did not match the actual distribution 
obtained by counting the number of the tubes from the networks. 
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One ofthe curves fitted to the network P curve is given 
in 6* 


irom tlae limited n'um'ber of analysis performed by ns, 
it is doubtful wlietber this model could at all relate the 


pore-size distributions to curves of real porous 

Irorn our study it can be said that networks 

usevl as a valid modol of porous media and pore— size 

could be obtained from P^ vs S curves. 

c w 


media, 
could be 
distributions 


*** 



CHAPTER 5 


CONCLUSIONS AHP RECOMvIEHDATIOH 

The classical model (bundle of tubes) of porous 
media litivo a number of weaknesses. As a result the inter- 
pretation of P^ curve based on this model is likely to be 
erroneous. Network is a better physical model of the porous 
media and P curves obtained from it is quite realistic. 

True poro-sise distributions (volumetric and number) can be bedAe-v 
estimated from P curves by modified Moyer's and modified 

V 

Eatt's method lv\<A\cx)A*A; dW wotk. 

Future VTOrk in this direction should ensure that the 

actual pore size distribution (number) of the simulated networks 

reasonably match the input distribution. More work is needed 

to find out the effect of tube length-radius ratio, the extent 

of interconnection, end number of tubes in the network on 

the P curve generated. Experimental v»rk with physical networks 
c 

arc necessary to verify these theoretical findings. 






Fig. la - A typical network used for generation of capillary pressure curves. 
V Input data table 1 . Numbers indicate tube rad. in microns. 

























Fig. 1c -Network from arbitrary bimodal distribution. Table 2. 











rig.le - Network used to obtain vs. 5^ curve ( Fig. 2e) indicating absense of the 
largest size, tubes in the input end. Input incomplete beta distribution 
with paraimeters 1*0,1*0r f^umbers indicate tube radius in micron. 





Fig. 1 f -Network for poisson distri 
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TABLE 2 

EBIQjRX EY DISTRIBUTIO N USEB IN SIMULATION OF 

NETWORKS 


ore-Size 

Cu'.ri. Number 

(,J,J 

Distribution 

0.0 

0.0 . 

2.0 

0.0434 

A.O 

0.1304 

6.0 

0.2608 

8.0 

0.3478 

10.0 

0.3912 

12.0 

0.4346 

14.0 

0.4780 

1 6. 0 

0.5650 

18.0 

0.9130 

20.0 

1.0 

0.0 ■ 

0.0 

_'3.0 

0.075 

,5.0 

0.325 

8.0 

0.400 

10.0 

0.450 

12.0 

0.500 

14.0 

0.550 

16.0 

0.650 

18.0 

0.950 . 

20.0 

1.0 


**** 
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table 4 


'’^1'ISON OF EXPj?RIMEH'I 

1)Il_.E0RE-SIZE DiOTRTBTTUTnAT 

11) ERO!"'; MERCURY PC 

iROSIliETRY TO TH4T 

COREECTEB 

BY MEYER' 

S METHOE 


C()ui't:s 'b ruujL , J 



I'o re -U Micron ) 

; :rt!M Ut i.pirm oqunl 

Lti 

Cum. Fractional, 
Pore Volume 
Exptl. 

Cum. Fractional 
Pore Volume 
Corrected 

1.066 

0.00578 

0.09797 

0. 6 3 '5 

0.01734 

0.03578 

0.213 

0.53757 

0.71379 

0. 1 066 

0.87861 

0.92197 

0.(}33 

0.93642 

0,93954 

■? — 

CM 

O 

0.99422 

0.99422 

0.01060 

1.00000 

1.00000 

Pyro: UP Fritted Glass 



1 'o rc - S i e ( M 1 c ro n ) 
Greater than equal 
to 

Cura. Fractional 
Pore Volume 

Exptl. 

Cum. Fractional 
Pore Vo lume 
Corrected 

1.066 

0.00592 

0.10029 

0.3d3 

0.07101 

0.24686 

0.213 

0.94083 

0.94083 

0.1066 

0.99408 

0.99408 

0.03331 

1.00000 

1.00000 

0.02132 

1.00000 

1.00000 . 

0,01060 

1.00000 

1.00000 
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'I'n.it'l « 4 (I'ontd), 

) -A-cti yo te d Clay s 


Po re - Si. ze ( Micron ) 
Greeter than equal 
to 

Cum. Fractional 
Pore Volume 
Exptl. 

Cum. Fractional 
Pore Volume 
Corrected 

1.066 

0.00513 

0.08693 

0. sot 

0.01026 

0.01604 

0.213 

0.4615 

0.14971 

0.1066 

0.36410 

©.44945 

0.05031 

0.63077 

0.63897 

0.02132 

0.85128 

0.85128 

0.01 060 

Itgjjyuy) opl;, 

1 . 00000 

1.00000 

J 0 r c; - : .1 i K c ( M i c r 0 n ) 
(;'i'(,'Mtc.' r Idrin Equal 
to 

Cum. Fractional 
Pore Volume 
Exptl. 

Cum. Fractional 
Pore Volume 

Co rrected 

1.066 

0.1685 

0. 16701 

0. 533 

0.22472 

0.51417 

■0.21:32 

0.60674 

0.75752 

0,1066 

0.74719 

0.79563 

0.05231 

0.85393 

0.85803 

0.02132 

0.94382 

0.94382 

0,01060 

1.00000 

1.G0000 
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TABLE 6 


C-ALCTJLATIOIVT OE TOBER Gilt. BISTRIBUTIOB MOM P_ DATA 
OE nET\?ORKS BY EATT’ S MODEL WITH MODIFIED EXPONENT 


Network 

P 


s 

Computed Number Distribution . 

G 


w 

V..= -2.5 

A..= -3.0 oi 

.=-3.5 

~~ 1V 


-Zl3) . 

“I4) 

(5) 

(Ti 

(7) 

I 

0.37 

14.00 

1.0 

1.0 

1.0 

1.0 


G.50 , 

10.5 

0.6920 

0.5247 

0.3854 

0.33 


0.75 

6.9 

0,5948 

0.3962 

0,26236 

0.201 


1.00 

• 5.16 

0.5155 

0.3086 

0.1883 

0.116 


1 .50 

3.46 

0.3685 

0.1716 

0.0904 

0.047 


2.0 

2.60 

0,2603 

0.0'31G 

0.011 

0.0044 

II 

0.324 

16.0 

1.0 : 

1 .0 

1.0 

1.0 


0.5 

10.5 

0.91 

0.831 

0.752 

0.667 


0.75 

6.9 

0.729 

0.558 

0.424 

0.313 


1.0 

5.16 

0.596 

0.390 

0.258 

0.161 


1.5 

3.46 

0.526 ' 

■ 0.315 

■ 0.195 

0.112 


2.0 

2.6 

0.318 

0,130 ■ 

0.064 

0.029 


3.0 

1.73 

0.235 : 

O'. 067 

0.027 

0.010 

III 

0.287 

' 18.00 

1.0 


V I' .O ' 

1 .0 


0.323 

16.00 

0.70 


0.485 

, 0.560 


0.370 

14.00 

0. 666 


. 0.43 

0.510 


0.430 

12.00 

0.666 


0.43 

0,510 


0.516 

10.00 

0,538 


0.287 

0.358 
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1‘able 6(contd) 


(2) 

(3) 

(4) 

(5) 

(6) 

0, 646 

8.00 

0.354 


0.1184 

0,862 

6.00 

0.283 


0.0678 

1 .290 

4.00 

0.221 


0.0365 

2.580 

2.00 

0.1015 


. 0.0 

0.259, 

20.00 

1.0 

,1.0 

1,0 

0.30 

17.4 

0.8657 

. 0.76 

. 0.72 

0.37 

14.0 

0.7195 

, 0.54 

. 0.44 

0.5 

10.5 

0.6296 

0.418 

0.305 

0.75 

6.9 

0.5143 

0.286 

. 0.180 

1.0 

5.16 

, 0,431 

0.204 

.0.118 

1.5 

3.46 

0.265 

0.07 

0.030 

2.0 

2,6 

0.260 

0.067 

0.029 

5.0 

1.73 

0.194 

0.027 

0.010 

0.37 

14.0 

1.0 

1 .0 

1.0 

0.50 

10.5 

0.5186 

0.215 

0.30 

0.75 

6.9 

0.5081 

: 0.202 

0.29 

1.00 

5.16 ' 

0.4032 

0.116 

0.182 

1.50 

3.46 

0.2763 

0.035 

0.07 

2.0 : 

2.6 

0.2206 

0.0128 

0.03 

3.0 

1.73 

0.190 

0.0043 ■ 

0.01 15 


( 1 ) 


(7) 


I? 


V 


0.16 

0.090 

0.042 

0.0 
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Table 6 (Contd.) 


( 2 ) 

( 3 ) 

( 4 ) 

( 5 ) 

( 6 ) 

0.259 

20.0 

1.0 

1.0 

1 . 0 

0.30 

16.0 

0.988 

0.98 

0.977 

0.37 

14.0 

0.8718 

0.795 

0.746 

0.50 

10.5 

0.6216 

0.4428 

0.360 

0.75 

6.9 

0.4009 

0.1818 

0.1144 

1 ,00 

5 . l '6 

0.3703 

0.1514 

0.091 

1.50 

3.46 

0.3029 

0.095 

0*054 

2.0 

2.6 

0.2077 

0.0275 

0.017 

3.0 

1.73 

0.1653 

0.002 

0*0056 



71 


HEPEHEHGES 

1. Eatt, I, Network, Model of Porous Media Parts I, II and 
III, Tran s.AIME 202 p. 144 (1956). 

2. Meyer, H.I., Pore Distribution in Porous Media, J. 

Appl. Hiys.,,2£, p 510 (1953). 

3. Dallavalle, J.M. , Microaeritics, Pitman, N.Y; (1948) . 

4. Carman, P.0,, The Plow of Gases in Porous Media, 
Buttorwortli Scientific Publications, London (f956). 

5. Dullien, E.A.L. and Batra, V.K., jDetermination of the 
Structure of Porous Media, I&EC, 6 ^ p. 25 (1970). 

6. Gates, J.I. and Teopelcaar Lietz, W, ; Drill, and Prod, 
Practice, API, N.Y (1950). 

7,, Eatt, I. and Dykstra, 'H. ; Trans. AIME 192 p 24-9 (1951) A 

8, Purcell, W. R, , Trans. AIME, 186 p. 39 (1949). 

9. Burdine, N.T., Gouroy, L.S. and Reichertz, P.P.; ^ 

Trans illME ISi p. 195 (1950). 

10. Gregg, S.J. and Sing, K.W, Adsorption, Surface Area 
and Porosity, p.139»Acadeniic Press (1967). 

11. Kozeny, J s S.B.ALkad.Wiss, Wien. Abt. Ila, 156 p.271 ( 1927). 

12. Carman, P.C, Trans. Inst, Chem. Eng. Bond. 1^., p.150 

(1937). _ , , 

15. Haughey,D.P. and Beveridge, G. S.G. ; Structural Properties 
of Packed Beds - A review,.- Can. J. Chem. Engg. 17), 

P.130 (1969). 



72 


14. 

15. 

16. 

17. 

18. 

19. 

20 . 

21 . 

22 . 

23. 

24. 

25. 


Kruyer,S. j 7r2ns.Pe.radaj Soc, 54 p« 1758 ( 1958) . 
iil’rtivel, L.J.V. find Iiressley, L, J. s iliial, Chem. ^ p.t492 
(1963). 

Meyer, _H.P, and Stowe, S.-A. r J.Colloid. Sci, 20, p. 893 
(1965). _ 

Moyer,' A P. and Stowe, H. A., J. Piiys.Chem, 70, p.3867 

( 1966 ).. 

Melrose, J.C. , V’cttability as Related to Capillary Action 
in Porous Media, SPEJ _5(3) 1965. 

Naar, J. and- Wy gal, R.J., Structure and Properties of 
Unconsolidated Aggregates. Can. J, Plays. 40 p.818 (1962). 

Markin, V.S,, Isv, Akad. Nauk.SSSR Ser, Khim. 1J,, P..1923 
(1963). _ y _ _ _ _ / - 

Dodd,C.G'. and Kiel, O.G. ; J. Pkys. Chem. ^ p. 1 646 (1959). 
Singhal, A. K., Partial Immobilization of Pluid Phases 
During Multiphase Plow Through Porous Media" ‘Ph.D. Thesis, 
Univ, California (Berkeley) 1971. 

Haring, R.E. , and Greenkorn, R.A. ; A:-. Statistical Model 
of a Porous Medium with Uon-uniforrii Pores" , A.I.Ch.E.J. 

16 (3) P.477 (1970). 

Pay takes, Tlon and Turiani A He?/ Model for Granular 
Porous Media, Part I,A.I.Ch,E. J. p. 58, Jan, 1973. . 

Ritter, H.L. and Drake, ,,„,Ii.C. , Pore Size Distribution in 
Porous Material, Ind.Engg.Chem. VT. (analytical edition) 
p. 782 ( 1945 ). 



73 


26. Zwie’fcefing, P, The -Structure and Properties of Porous 
Solids, p.287 Butter?/orth Publication ( 1958 ). 

27. Dullien, i'.A.L,, Determination of Pore*Accessibilities 
- an Approach, J.P.T. p. 14-15, Jan. 1969 . 

28. Levorett, M.C.; Capillary Behaviour in Porous Solids, 
Trans.iLEME 142 p. 151 (I9|-l).. 

29 . Leverott, M.C., Lewis, 17. B, and True, M.E,, Dimensional 
Model Studios of Oil Field Behaviour, Trans.AME 146 , 
p. 175 ( 1942 ). 

30 . Hose,W. and Bruce, W, A,; Evaluation of Capillary 
Character in Petroleum Reservoir Rock, Trans.AIME 186 , 

p. 127 ( 1949 ). 

31 . Fatt, I., In Discussion of Paper, 'Evaluation of Monte- 
Carlo Methods in Studying Eluid-Fluid Displacement and 
\7ottability in Porous Rocks* by C.G. Dodd and O.G. Kiel 
(Ref, 21 ) , 

32 . Wienbrandt, R.M. and' Fatt, I. 1 A Scanning Electron _ 
Microscope Study of the Pore Structure of Sandstones, 
JPI, p.545 May 1969.' 


*** 



74 


APPMDIX A 

NETWORK MODEL s PLOW COMPUgAPIOl 


Pressure Oomputation s 

Assuming that flow in tlie pore cliannels of a porous 
system is essentially viseous in nature (Re 4 75), Poiseuillels 
law of viscous flow could be arplied for flow computation in 
each tube of the network. For single phase flow in a cylindrical 
tube of length ond radius rj_, the expression for flow rate, 
is given ass 


or 


~TJZ 


AP 

h 


(A.1) 


4i, = 

v^here 4SPj_ is the pressure differential causing flow in the 
ith tube. For multiphase flow and h.,Pj_ may be modified 
depending upon the relative proportion of the different fluids 
present and the nature of the flow regime respectively. Knowing 
total fluid conductivity, K^, of each tube and maintaining total 
volumetric balance a Kirchoff type equation can be written for 


each node point of the network; 

&I< h,3 - A,i) + = o 

k=i“1,i+1 

(i,0)5^k!ir^ 


(A.2) 



75 


Knowing the values of K^, and the pressure at the 
hounda-ries, the set of sjjnultaneous algebraic equations can 

be solved by Sauss-elimination technique. 

Ganillarv Pressure ; 

Capillary pressure, P^, is applj.cable only to those pores 
which undergo displacetnent flow. Por .single phase, slug or 
annular flov/ regimes, caiiill ary pressure is essentially zero. 

In this computation P^^ is taken as equal to the pressure required 
to squeeze an fluid/fluid interface from a mode point op>ening 
(radius r^^) to a tube of radius r^(rjjQ t ^) » Por a tube of 
triangular cross-section 

= 1.78 Y^ooa 9 (-! a) (A.3) 

^no n 

where = interfacial tension and © is the contact angle which 
in our case is the receding contact angle. 

Plow Behaviour in a Single Tube s 
i) Single Phase Plow; 

Iho total volumetric flow rate through a tube of 
triangular (equilateral) cross section, with the radius of the 
inscribed circle, r, viscosity of the fluid /:x, and pressure 
gradient iiP//j>l, is 

g _ liJ^ , eI. ( ^•4-) 

: jr 'A ^ - ■■ 

ii) Annular Plow; 

We used modified luster relationship as suggested 
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22 ’ 
Toy Singhal for flow computation. 

,4 


q. 


w 


% = 


1.35 
• v/3 
1.35 


r 


i.\P cP 

— % 


ZM 


tft 


M. 


n 


2S + 


(A.5> 


^{1-2 


A 


£))1. 

Kir 

w / . 


(A.6) 


whBTQ S d8Bo tos Batur^tioB and subscirip^s w and n denote Vv^etting and 
and. non-wetting phcose respectively. 


(iii) Slug Plow; 

Total flow rate is given by 


q 


1.33 


tot rf 4^ + '^n^n + 


ziP 


V'' VV iJ U w 

where a is the side of the triangular cross section. 
( iv) Displacement Plow: 




Selection of Tube for Plow Ohange ; 

Depending on the availability of fluids at a node point, 
transition between flow regimes will occur. In many instances, 
channels leaving a node point may be capable of carrying more 
of a particular phase (fluid) than is available at that time. 

In such cases, saturation and flow regime changes in the tubes 
will have to occur to unload all available fluids. A table 
of possible transition between flow regimes is given in Table A-1 . 
Some times there may exist more than one pore which might be 
undergoing flow in the regime to^ be altered. Selection of tubes 



77 


for flow changes in such cases is done by Monte-Carlo method 
with a bias for flo?; capacities of the competing tubes. 

General Scheme for flow OocnDutation ; 

Yfortis in block capital letters are subroutine name 
referred to in the program listing. 

The network is generated in the program by SIMU with 
specified rx)re-size distribution and grid size. This subroutine 
also computes the volume of the pores connected to each node 
point and flow conductivities to the two fluids. ■ Other data^ 
fluid viscosities, interfacial tension, pressure differential 
across the model and contact angle, are read in the main^^^ ^ ^ ^ ^ , 

m “ . . , ■ 

program CALC, , The radius of the opening at the node points, r^, 

is chosen to be slightly larger than the largest tu'Se radius 
of the network. 

The next step is to determine the pressure distribution 
in the model. Pressure at each node point is computed by 
solving the system of equation A-2 by PRSOI through SOLYE. 

Apart from updating the conductivities of each tube after every 
time step, SOLYE also immobilizes those tubes which may cause 
outflow from the first row and inflow from the last row. Also, 
those tubes, whose eapillary pressure is greater than the* 
pressure differential available and are oppositely directed, 
are blocked by the action of trapping factor. Y/hen a tube is 
blocked its conductivity is set equal to zero and the system 
of equations A~2 is solved* This is repeated till no more 



78 


tubes are blocked. While calculating T for tubes in the 
first row, the node point opening is assumed to be infinite. 

I'or other rows the specified value of r^ is used for the 
calculation. 

The total outflow' of each phase (fluid) from each node 
point is calculated under the existing pressure distribution 
by PLOVi/, Outflow from a tube is entered as inflow to the 
connected node point by 00^(1). The imbalance between inflow/ 
and outflow of a particular phase at a node point is corrected 
by CHMGE which also effects the flow changes necessary for 
this correction. The capillary pressure and flow conductivities 
for the corrected flow distribution is computed by C0ED(2) and 
stored.. Outflow from the model is obtained by adding together 
the inflow to the last row of node points. The volume of the 
fluid left in the model is computed to evaluate the saturation. 
This process is repeated till the equilibrium wetting phase 
saturation in model is obtained for that pressure drop. 

Steps of Computation ; 

(1) Call SIMU to simulate the network model. Assign 
tube sizes, locations etc. Calculate the volume of each tube, 

. its conductivity to displacing and displaced fluids. i . . 

(ii) Call SOLVE to update the conductivities through 
G0ND(2) and solve the system of "equations A. 2 hy PRSOL to 
evaluate the absolute pressure at each node point. In case 
some of the tubes are blockedf 'call C0NI)(1) to zero the 
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conductivity of those tube s and again call PRSOL to solve 
the systemof equation A-2, Repeat till no more tube s are 
blocked. 

(iii) Choose node points in descending order of absolute 
pressure and call PLOW to calculate outflow of each fluid from' 
each tube and node point. 

(iv) Calculate imbalance at each node point to test 
whether overall material balance of incoming and outgoing fluids 
is satisfied. 

(v) Imbalance in the material balance of each phase is 
* removed by assigning the excess incoming phase to various 

tubes by Monte Carlo method through CHMGE without forcing flow 
change in the tubes. 

(vi) Assign- outflow from each tube as the inflow- to the . 

connBcijOd nod© points# 

(vii) In case of back flow, calculate inflow to those 
tubes and node points in PLOW. If there is any imbalance either 
in total fluid at the connected node points or imbalance in a^^ 
phase in any tube, call GHMGE to assign the differential volume 
to various tubes. 

(viii) Calculate outflow from the model by adding the_ 
inflow to each node points in the last row of the network. 

(ix) Calculate total fluid volume left in the model by 

n j>- 4«vvi<~v -PI n *1 ini Q-i "tli0 0nd oi "til© 

adding the volume qf the fluid in eacn uupe 
time interval.' ^ ■„ 
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(ac) Calculate the saturatioh for each phase. 

(xi) Repeat steps (i) to (x) till there is practically 
no change in the total outflow from the model in the successive 
time steps (less than of the outflow in the first time step) 
for the same ^ P. 

(xii) Test whether the sum of the cumulative outflow and 
the volume of fluid remaining in the network for .each fluid 
at the end of every time step is constant. 

(xiii) Repeat the sequence of operations with higher AP, 
to get the value of eguilihrium (wetting phase saturation) for 
different Pq. Thus curve is generated. 
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TABLE A-1 

TILmSITION BET?.rEEM ETiOW BErTTMTj'.F! 


Plow Regimes 

To / Prom 

Single 

Phase 

Displacement 

Slug 

Annular 

Single Phase 

X 

1 * 

1 * 

1 

DisPlac©tnent 

2 

X 

- 

- 

Slug 

- 

- 

X 

1 

Annular 


1 *-^ 

1 

1 


X No change in flow regimes. 

1 Easy transition* 

2 Transition influenced by capillary forces, 

.Transition not possible without going through 
other flow regimes, 

* Modified single phase flow for non-wetting flow 
because of a film of wetting fluid at the pore 
walls, ■ 

** Possible only' for non-wetting phase displacing the 
wetting phase. 





82 


APPENDIX B 
METER’ S MEIHQD 


The radius r of a pore is defined to be equal to the 
radius of a capillary tube which will just alio w mercury to 
enter it against capillary force at the same pressure as that 
at which the pore v/ill fill, if connected to a mercury source. 

Lot be the fraction (volumetric) of the total 

void apace occupied by r-pores-i.e. pores whose radii lie in 
the range r to r+dr. The. number of such r pores in a unit 
volume will be given by ■ 




if \ ■ 

a(r) dr = — . * (B,1) 

Kr^' 

where f = fractional porosity 

Kr - vol. of a single pore of radius r 
It is necessary to designate one representative point 
in the interior of each pore which we shall term "centre of 
the pore". For our purpose, any point will do so long as it 
is not too close to the wall of the pore and is consistently 
chosen, o.g. c^g, of the void space will serve. Also, we define 
a volume dv so small that the probability of the centres of 


two or more r-pores being located in this volume may be 
neglected. 

Let the probability of ra r-pores being located 
to any volume v be distributed pores 
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dr = -(adr) dvl 

'^-^_1 (v,r)dr (adr)dv 


(B.2) 


wliich on re-arranging gives 


dr) r(a dr)vj - - 

— _ „ exp ^-(adr) vl (B*3) 

Now we define a volume v(rQ,r) such, that, if an i-. pore 
is located in v(rQ,r), it may be expected to be connected to 
an rQ-pore. This can be done if we assume that two pores are 
connected if the centre of the smaller lies in a volume 
adjacent to the larger pore, which is the difference between 
the volume, of a pore whose radiu#'-‘equals the sum of the radii 
of the two pores minus the volume of the larger pore i.e^ for 

^o> ^ ' 0 

v(rQ,r) = K(rQ+r)^ - Kr^^ = K(3r^ r^ + 5r r^ + r^) {B.4) 

Now if in Eq. ( B-3 ) v is replaced by v( r^ ^ r) , we get the\ 
probability of ra r-pores all being joined to the same rj^-pore, 
The probable number of such groups of pores which would exist 
in a unit volume would be given by the product of. this pro- 
bability times the number of v (rj^,r) volumes in a unit cube 
(which is oqu.cl to the number of r^ pores in a unit cube). 
This number is equal to , ■ ^ 


a(r- )dr .(wtej. 


■-‘fj 


a(rQ)dr 
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We are here interested in r)dr which is- the 

probability of an r^ pore not being connected to an r- pore. 

Let a porous rock be subjected to mercury pressure such 
that all pores of radius ^ o? would fill if they had access to 
mercury, and let f (r) be the measured fraction of the pore 
volume that is filled up with 'injected mercury at that pressure. 
Lot 0(r) be the actual fraction of the pore-volume devoted to 
pores of radius ^,r. 

If the simplifying assumption is made that only pores 


of radius r which will fail to fill are those that are not 


connected to any other pores of equal or larger radius* That 

is the probability of an r pore failing to fill is given by 

the multiplicative lav/ of prl5bal>|j|^'ty (r ,r)dr. 

•<' r . ■ ° 

The total number of pores of radiusj^ rg which will fail 
to fill is given by the expression 



Hence \ (r) 


= r a(r)dr 

sj' 

■ ■ ■ ' • , 

(■B.5) 


T , ■ ■■ 



'• 1'' ( r. , r ) 1 


= 0(r) [1 - 

i 0 i 

. J ■ 

(B.6) 

J a(r)dr * 



This is the expression from which we shall obtain true pore- 
size distribution ©(r) from the actual distribution f(r). 


Numerical Approximation ? 

Since the above expression is not convenient for numerical 
calculation the following approximation are used. The radius x 
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is disci'i'tizsd so ilial; it takss only a finite nmalDetr of 

values and 

Z.\r 

= 

r. 


^i-1 “ ^i 

j'”'^ a(r)dr 


1,2 


“ i (r^) \ ) 

'tlj = ®^p l-=»i 

,t(ri,rj)= “k 7T 'Pkl 

k =1 1=1 


% 


Substituting these expressions inEq.(B. 6) gives the approxt- 


, ' d £ fr f ki 

. JL.. ^ * r. 1^1 ^ 1=1 1 Vp f:„’) 

k = '^JA ^ U, - “r" i ‘ 

' ■ ■ 

■ k=r 


mate equation 

n iXi A* ^ (to** *1^ ri 

r/) = 'T v'l. = ^C"-)y,_ X i1 - 

k=1 k=1 


Algorithm ; 

'( 1) Assume a value of and calculate and 
■ (2) Using these values solve Eq. (B. 6a) for l/^, 

(3) If computed J;. value is not close enough to the 
assuni6d'value, mkke another assumption and repeat 
this proc'^s. 
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APIMDIX C 

GLOSSARY OF TEMS USED 

Porosity (f) is the ratio "bebweGn the void srace and bulk 
volume in a porous medium. 

Q.a ilil 1 ary Pr o s su ro (P^) .is tho difference in pressures across 
a fluid/fluid interface. Capillary pressure curve describes 
the changes in capillary (threshold) pressure of a pocrous 
as a function of its fluid saturation. 

Lhre sho Id Pro s sure (L^) is the minimum excess pressure must 
be maintained in the non-wetting phase to enable it enter 
a given porous medium. 

Satu ration ( 8 ) of any fluid in a porous medium is the ratio 
between the volume of the ‘'fluid in the medium to its pore 
volume. 

Contact ifigle (Q) is the angle (through the displacing phase) 
th.ot the fluid/fluid interface makes with the surface of the 
porous medium. 

Wotting Phase is the fluid in a porous medium ?*ose contact 
angle with tho solid surface is less than 90 '^. If the contact 
angle is more than 90® it is called non-wetting phase. 


**** 
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N 

SIGMA 

THEIA 

POR 

R(I) 

PC 

SGAMA 

BSTKU 

NU 

DIMU 
SV(I,J) 
PSI(I, J) 
PSIR(I, J ' ■ 


APPEIDIX D 

PROGBAI'I LISTING OP mJER' S MEIHOP 

lumber of data in the set 

Interfacial tension 

Contact angle of non-wetting phase 

Porosity 

Pore-size corresponding to P(I) 

Capillary pressure 

Cumulative Volumeof Mercury Injected. 

Assumed value of nu 
Computed nn 

Difference between estimated and computed nu 
v(r-,r-) s characteristic volume 



ISN 


scuRce STAT£f^e:a 


FORTRAN SOURCe LIST 


m 


0 SIBFTC f^AIN 

C MEYERS METHOD 

G PC IN PSI 

C CAT*1. MEANS THE USUAL VALUES OF CHANGE ARE NOT TO BE USED 

1 COMMON/Al/PSI(20,20I,TEMP(20,20),SUMA{20)iSGAMA{2C} 

2 CCMMON/A2/RS{20)tRC(2 ), SIGMA 

3 C0MM0N/A3/NU{50’'j) , FOR 

4 DIMENSION P(20) tSV ( 20 ,20 ) , R ( 20 ) , A { 20 ) ,D IFNU ( 500 ) , E STNU { SCO ) 

5 real NU 

6 READ 6u ,N . : . 

lb 6u FORMAT (14) 

11 PRINT120 

12 120 FORMAT (IHl ) 

13 REAC62, SIGMA, THETA 

14 62 FORMAT (5F10. 3) 

15 65 FORMAT (5F10. 5) 

16 THETA® TH£TA*3. 1416/180. 

17 REAC62 ,PQR 

20 PRINT50, THETA 

21 50 FGRMATiSXj’kTHETA®*, FIO.5} 

22 N=N+l 

23 0012 I =2,N 

24- R(I ■ ■■ ^ . :: 

25 RS ( n*o.o 

26 RC( I)=0.0 

27 A(n»Q.O 

30 SUMA (I )=0.0 . 

31 DO 12 J»1»N ' , 

32 SV(I,J)"0.*j ■ ■■v..,;;,,, . 

33 PSI(I,J)=C.O 

34 PSI R =0.0 

35 12 CONTINUE ■: 

40 00 101 IK=2,N 

41 101 REAC65 ,P( IK) ,SGAMA( IK) 

43 PRINT20, (PdK), IK=2 ,N), (SGAMAdK), IK=2,N) 

54 DC 40 IM=2,N 

55 4C SGAMA( IM)®SGAMA{ IM)/PCR 


57 

60 


P R I fs T5 

5 FORMAT ( 5X, 12HACTUAL VALUE ,10X,2HIK 


,5X,10HTRUE VALUE 




OOiC J»2,N . ' . 

CAT»G.O ■ . * 

CIFNU(l)»b*0 ^ ^ 

I K® 2 • f ' 

ESTNU( IK)»SGAMA( J) i , .t 

2 GALL S n I K , P S IRt A # ! t ? ’ f ? L i 

IF(IK.eQ.2.AND*’s> . f fRI NT20 » Jfi 1 1 1 » N ) 

I F ( I K. eo. 2 ♦ AND* I J > , 1 0 





CGG145 


ISN 


133 


136 


141 


144 


147 


152 


155 


156 


161 


162 

25 

163 


i,66 


171 


172 


173 


174 

7 

175 

6 

176 

i«' 

200 


205 

8 

206 

22 

207 



.CGGl 4 5 


5CUPCE STAT^^'ENT 


FORTRAN SOURCE LIST MAIN 


I F( ESTNU( IK) .GT .O.OCi )CHANGE=C.001 
IFIESTNUdK) ,gT.' .• = l)CHANGE=fl.Ci 
IFtESTNUdK) .GT.v.l)CFANGE=G.l 
IF(ESTNiJ( IK) .GT.0,9)ChAIM6e=0.C05 

1 F (C I FNU{ IK ) ,GT .0.0 ) ESTNUl IK+ll =t STNU { IK ) + CHANGl 
I F (C IFNUI I K) .LT.f ,!’ ) LSTNU ( I K+1 ) = E STNU { IK ) -CHANGE 
IK = IK+1 

IF ( IK.GE.5C0) GO TO 7 
GO TO 2 

CHANGE = CHANGf:/2 . 

I F ( Cl FOiJ ( IK) .GT .0.0 ) cSTNUI IK + 1 )=eSTNU{ IK )+CHAN6E 
I F{ CIFNUI IK) .LT .0,0 ) tSTNUI IK+1 )=E STNU ( IK )-CHAN6E 
IK* IK+1 
CAT»1.'J 
GC TO 2 

PRINT 6iSGAMA(J),IK,eSTNU(IK) 

F0RNAT(5X, FlG.5,ltX, I5,5X,F1C.5) 

CONTINUE 

PRINTSf {A( J) ,SUMA( J),J=2,N) 

FGRRATdOX.ElS, 7,iOX,Fl5.7) 

STOP 

LNC 

IBMAP ASSEMBLY MAIN 


NO MESSAGES FOR ABOVE ASSEMBLY 


« 


CGGi45 

ISN 


SOLPCt STA^r.^ENT 


FORTRAN SOURCE LIST 


C' 

1 

2 

3 

4 

5 

6 
7 

10 

^ 11 
12 
13 

15 

16 
17 
213 
21 
22 

■ 25 
26 
27 
30 

33 

34 

35 

40 

41 

C6G145 


$IBFTC SL0I 

SUBRUUTINl: VOL ( P ,R , IHLTA, IK, N , SV ,M , IM } 

■ CONPriN/Al/PSI (20,20 ),TENP( 20,20 ),SUMA{ 20} fSGAMA (20) 
C0MNri!vl/&2/RS{2r ),rC(2 .),SIGMA 
CCNNQN/A3/NU(50.'; ) ,FCR 

DIREMSIUN P( 20) ,SV( 20,20) ,R{20) , A (20 ) , OIFNU ( 500 ) 
FACT0R=1 ./6U947.3*1C.««4 
on 1 I=2,K 

R( I )=FACTCR’!'2.*SIGMA»!'C0S{THETA)/P( I ) 

R( n=ABS(R(I ) ) 

RS( I ) = R( I )*’!'2 
1 RC(I)=‘R(I )=^=!'3 
R ( 1 ) *2'. • A 
RS (1 )«R( I )442 
N1»N+1 ^ 

DC « I»2,N 
CO 6 J*?.,N 

I F (R ( J) ,GT.R( I) )G0 10 7 

5 SV( I,J )»2.*R(n *R( J1+RS(JI 
GO TO 6 

7 SV( 1,J)=2.4R(I)*R( J)+RS(I) 

6 COMINUB 
DOIC I=2,N 
DClC J*2,N 

IT SV( I,J)=A8S(SV( I,J) } 

RETIRN 

END 

IBMAP ASSEMBLY SUBl 


*=10 


NO MESSAGES FOR ABOVE ASSEMBLY 


^1 


CGG145 

ISN 


SCURCt STAT'^fiENT 


rj $IBFTC SUE.: 


FCRTRAN SOURCE LIST 



assembly SUB2 

.‘■i' ■' ■ 


CGG145 


" NO MESSAGES FOR ABOVE ’ ASSE^to, 
C«G145 ' v> 4,'.,,'’ 


OO.pl)Oli 




*** ' OBJECT 


SLBPUIITlMi- Sl{lKtPSIR,A,ESTNU,THETA,N,P,R,SVf IL,JKr 
CUMf'flN/Ai/PSI (2 t2sU,TEMPl20,20},SUMA(20),SGAMA{20) 
COKMON/A2/HS(20),RC 120), SIGMA 
CCMMnM/A3/KU(5'-V’ ) , fCR 

Ci MENS I UN P( ?*)) ,SV( 2O,20),R(2CI,A( 20 ) , DI FNU ( 5CC ) ,ESTNU{5( 
REAL NU 

CALL FLUN{G9S99) 

N 1 =N 

NUUK) =ESTNU{IK ) 

SUMA (1 ) = 0.U 

IF( JK.EQ.2)CALLV0L ( F,R,THETA, IK,N, SV,M,IM) 

I = IL 

A( I )«NU( IK )*(1./R(I )-I./R(I-13)*P0R 
DC 12 J=Z,N1 
B=-A(I )«SV(I , J) 

GO TO I 
GO TO 2 


IF(E.GE. 88.03 
IFiE.LF. (-S8.0) 3 
PSI (I, J3«E)(P{B) 

GO TO 12 

1 CONTINUE 

3 PSI (I,J) = l.uC30 , 

GO TO 12 

2 CONTINUE 

4 PSI (I, J3aO.U V 
12 CONTINUE 

J = JK 

DO 8 K*2,1L 
T EM P ( I » J 3 -i . ' 

pRCC*i. ' ■ V 

DO 10 L=2,JK , ■■.■'■-..■A- . 

0ELTAR=R(L-1 )-R(L3 
PRCC»PROD*PSI (K,L) 
IF(PROO.LE.O.ie-373 FR0O=0.0 
10 CONTINUE « 

8 TEMP (I ,J) “TEMPI I , J 3 +A ( K) *PROD 

PSIP“TEMP{ I, J) tv 

9 CONTINUE 

SUMAd )»SUMA(I-1)+A(| 3 , 

5 F0RMAT(5X,8E15.73 


RETURN 

END 


' 21 MTS. 
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Avmwix E 

PR0GEiU4 LISTING gETWQRE MOD'RT, 

The first subscript in each variable in the program 
refers to the node point. 

X(l) to .X(8) The radii of the 8 tubes 
X(9) to X( 16 ) Length, of the tubes 

X(l7)to X(24) Fluid content and flow regimes in the tubes 

1.0 wetting (one phase flow) 

■ t,. 

2.0 non-wetting (one phase flow) 

3.0 Displacement, (non-wetting close to the 

node point) ■ ’ 

3.5 Displacement (wetting closer to the node 
point) 

4.1 Annular flow (displaced phase wetting) 

4.2 Annular flow (displacing phase wetting) 

4.5 Slug flow 

X(25)to X(32) Flow or no flo?/ situation 

1.0 Flow occurring 

2.0 Flow blocked 

1.5 Flow provisionally blocked 

2.5 Flow permanently blocked 

X(33)toX(40) Flow conductivities of the tubes to displaced fluid 

X(4l)toX(48) Flow conductiTtties, of the tubes to displacing 

fluid 

X(49)toX(56) Yo lume O'f the "tubes 

X(57)toX(64) Displaced fluid '^^teut of the tube 
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X(65) to X(72) 
X(73) 


X(74) 

X(75) to X(82) 
X(83) 

X(84) 

X(85) to X(92) 
X(93) 

X(94) 

X(95) to X(102) 

X(l03)to X(110) 

X(111) to X(118) 

X(119) to X(126) 

IB 
JR 
GMA 
THETA 
RAD 




Effective pressure drops across the tubes 
Incoming fluids 

1.0 displaced 

2.0 displacing 

5.0 Both 

4.0 none (node point isolated) 

Absolute pressure at node points 

Appropriate conductivities of the tubes 

Volume of incoraing displaced fluid 

Volume of incoming displacing fluid 

Effective capillary pressures 

Provisional total outflow of displaced 
fluid 

Provisional total outflow of displacing 
fluid 

Provisional inflow of displaced fluid in 
tubes 

Provisional inflow of displacing fluid 

in tubes 

Provisional outflow of displaced fluid 
from tubes 

Provisional outflow of displacing fluid., 
from tubes 

No. of rows,, in network 
No. of coluton in network 
Interfacial tension 
: Contact, angle 
'>•' Tu'lje' rod.ihf.; ' 
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COM 

Cumulative distribution 


GRID 

Grid size 


TRAP 

Value of trapping factor 


HYST 

Hysteresis in contact angle 

(=0) 

HI®; 

Time step 


ISTAR-IFIN 

Ho. of time steps for each 

pressure 

AV 

Displaced fluid content of 

the net?;ork 

ICP 

No, of pressure steps 




STATL'n xjT 


SOUR CL 


FORTRAN SOURCE LIST 




$IBFTC calc 3PP 

'CLiHMON/'-Yl'InKC,/'' ( F , , If ) , AH(5''' ,8) t GC , VO, VW , TI YL , VJ F < 3 ' , 3 ) , I J , I K , JR , 3 P P 

j.XA ,!■("( ) ,HYST,CAT,3AMMA, JjFLOtS: ,8} ,GR1D 

CUHI'1CM/''>K*./'S/I1M, A( 5 0,30) , AB( 50) SPPOO'.j' 

COMM .IN/SOL/a/' { 3 5 ' ), dPI 5 j ) , AVI 4D) ,WR( A n , W0( 4. ) , W1 { A ) 3?Pp 3 

COMM 1M/'..JL;VC AP {'■. i, *,-), 1, I CHANG, ARAT, TRAP liPP,. 7 

CUMM;)N/C ilA''./i.!i F j, WO IF .^P‘'>00.), ' 

CUMMOJ/CHi J/ICA vi(50,B),AWDIF{50,8),A01F0(50,C.),BASL 3 p^^'J j.ri'.. 

CUMMr)Nl/3'JLC/ALj3,CP, JIM, ICP 3PP„._ .. 

CQfiMDM/Sl M/C JH( ) 5) , RA0(x5) ,KL ,LK HpPC ! 

C Rt AD AND PRINT ALL TO"- RELEVANT INPUT DATA 3 PPO' 3'.; 

CALL FL0N( 09'0)9) RpP/ .. L , 


DO 9 .. lKj=j, 

read ! , IR,JK, GAMMA, THETA, CP, VO, VW 
PR1NT9 yi 
93/ FOKMAUIH;. ) 

PRlNl ^ , Ik, JR, gamma, theta, CP, VO, VW 

READ ;U,KL,LK. 


3 ppr x;. 
dPP*' ' 
3 PPOjx 3.. 
3PPM 13., 
RPPODi A'- j 
3 ppOOlbu 


PRIN'i iJ, Ii:, JR, KL,L K 
read , (CIJM( I ) , I = I ,LK) 

PRINT 9' , ( COMd ) , I = I,LK) 

READ 90, (KADI J) ,J = l ,LK) 

PRINT 9 , (KAU( J ) , J=i,LK) 

PRINT 9 ,ViJ,VW,CP 

READ 9',),GRiD 

PRINT 90 , GRID 

READ jL2,XA,HYST,TRAP 

PRINT 12, XA,HYST,TRAP 

READ 7, TIME 

PRINT 7, TIME 

READ I..,ISTAR,IFIN 

PRINT 10,ISTAR, IFIN 

I =I STA R-i 

READ 37,WR(l + i) ,D£N (I+1),STW,ST0 
IF( ISTAR.EO.l ) 1=1 
READ 37,AV(! ),WO(i) , WOID.WKI) 
A00=W3( I ) 


7 FORMAT! iX,Fl5.1o) * 

10 FORMAT! 213, 6FiO. 5) 

11 FORMAT! 20X,4I9) 

12 FORMAT !8FIC. 5) 

13 FORMAT! 6FI..U 5) 

14 F0RMATI12F10.5) 

40 FORMAT! IX,0F1G.5) 

3,0' FORMAT! IX, 21 8) 

50 FORMAT !20X, I5H0I L S ATURAT I0N=,F13. 4, 17HWATER SATURATI 0N=, FI 3. 

51 FORMAT! 20X,7HTH£ PT.,I4,25HIS BETWEEN lAND 4 IN C73 > 

52 FORMAT(20X,7HTHe PT I4,12HMAS 4 INC73; 1 

100 FORMAT!//// 8!9X,iH )l "C ^ . 

Ul F0RMAT(30X,21HPRESSURE OISTR^.S^TION T 

117 ‘FORMAT, (//IX, 3uX,l8HFLOtO., il$T,R:i8UTI0!*l ^ 

759, „EmMAT C20X, ^CONTACT - 

7Sa '.format (2CJX,*VISC0'SI,':T,Y' 


b P P 0 ■ ' 1 L 

BPP'J x7 

B P o 00 1 o j 

B PPOOI 9C 

3PP0‘'? ; . 

6PP0 ~2i ■ 

BPP0022G 

3PPrj:;'23' 

3PP0/29.. 

BPP0':»2 5 

3PP0026G 

BPPQu27c 

BPP00280 

BPP00290 

3PP0O3'Jl 

8PP033i-^ 

3PP00320 

3PP0033J 

BPP0G34^ 

BPPOGaS'* 

3PP00360 

3PP003?i 

8 P POD 380 

BPP00390 

3PP00400 

9PP009i;> 

SPP00420 

3PPCr'430 

BPP0:‘44.; 

3PP0045D 

BPP0O46D 

BPPn^'-A? ' 

DPP'T.^b. 

BPPiJ "(99 

aPPODBO. 

app<; , 5 

3 PPG '57 


SO'JRCI-: RT/W - .v'i _xn 


fortran source list calc 



. .J 


FORM'-' T ('ll.) 

Jli'‘.=0- .JR 

THI;:TA is CiiM'Ow'i AnIGL.' AND GAMMA IS IMTERFACIAL TENSION 
ARAT=V'0/VN 

ARG=TH r T A* " . 4?, o/ 1 8C« 

GC = 1.7o.’*GAMnAN'CI}S ( ARG)/bS947.31 
CAT=C. D 

PRINT 759, TH 'O'A 
read data from TAPE 

call si mu 

PRINTA7,( CM J,<),K=l,d),J = i,N} 

PRINT3 7, ( (XM J,<) ,K=57,64},J = 1,N) 

DO 2 ’ lU = i f N 
1K=IJ/JK+.. 

JK=I IK-^ 1 
1F( JK.Lr.'O) K = IK-i 
IF( JK.Llf . MJK = JR 
EP( IK, UK) =)M I.J, 74) 

CAP(IKtUK)= . .1 
DU 39 IA=:. ,'i 
X(IJ,IA+16)= .3 ■ 

X( lJ,IA+24)= ‘ .9 
Xd J,lA + i36)*Hl*TA 
1F(X(I JjIAj.fQ. lJ,IA + 24) = 2.j 

X( IJ, I A + 40 )»X{ I J, IA+32)=^'(/0/VW 
X(I J»I A+56)=X( I J,iA+43) 

X{IJ,IA+74)=X(I J,lA+32) 

X(IJ,I A+84)= .3 
X( IJ,IA + 94) = n.3 
X(I J,I A+102)=0.0 
X{IJ,1A+1I )=.>.„> 

X(IJ,IA+llU)«u*.U 
39 CONTINJi: 

XI I J, 3 3) =0.3 
X{ I J,84 ) =y./ 

X(IJ,93)=0,C 
X(IJ,7-i)=1.0 
X( 1J,94)=0.':'‘ 

2D CONTINUE 

PRINT 12, GAMMA 
STW = ('i.fJ 

ST0*3. 3 ; 

AV(l)*U.u 

."'/.DO 900 ICP = J. ,12 
M'Wn l)«U.O 
DO 903 ICP»I,I3 
IF{ ICP.SE.'i)TIMe»O.0Ol ' . ' A 
IFnCP.GE.5)TIME»0.005 ' 

IF{1CP.Ge.8)TIME=O.Ol " ‘’C'.,.. ,v. . 

READ7,CP ' '-‘M Jit •*/'':'■ '''x'.-'t 

PRINT?,' CP ' :t. 

DO "too I«rSTAR, IFIM ’ 

I^'CICP.E'Q. 3 ) 

'iFticP.EQ. i) wi { I , "■ , 


XPP ’ 5 

3PP’’ 

5 PPO tSi A 
3 H P t 
APPV 57 , 
ePP-.i P'.j 

3PP-; A . 
‘"ppu e: 
SPP0'')6..U 
3 PPO'Cb '.'O 
b P P L . f> - j 
A I-’ : ' ;j> 3 2 

3 PP0-!6‘^.:.i 
iPP 00.55 0 
3 P P \i o 6 ■ 
BPP'; D? I 

appooooo 

3 PPOO^sUO 
3PP0 .7 
BPP A -7 i . 
appol) 
apprj' 7 5 .■ 
BPPO. 74.- 
8 PPO,. 73 J 
BPP00760 
BPPQ.-77.- 
BPP0>'7B,< 
BPP9C790 
Bppooaou 

3PPCI0810 

BPPQOSkA) 

BPPOQ830 

BPPOOB40 

3 PP 0085 a 

BPP0T36D 

BPP00870 

Bpponefao 

3PPOO&9’i 

BPP0>:t9U',.» 

BPP00920 

3PP00930 

3PPO094U 

dPP0395U 

3 PP 00 960 

8 PP 00 971' 

3PPCiri9ao 

3PP0U98i 

BPP0U98i 

BPP009fe5 

BPP0'..99. 

bPPO;.OOo 

BFPOJ Ci u 

BPP’ 

R P P 0’ *, ‘ 3 


^7 


FORTRAN SOURCt LIST CALC 


SOURCE STATE Mt. NT 
41a CQUTtMUE 

: FRaM - vT'AMlEj': A.^D SOL Vc THE SYSTEM QF 

CALL SiUVc 

1 F { ICP .LT. C ) SO TO .. ,>7 
).'!8 PRINT : 

: PRINT PRCcjS DiSTRIBLTlQN 

DO i07 IB=. , IR 
PRINT . 

PRINT .Tf ( . P ( n , JK ) ,JK=i, JR) 

IJl CONTINUE 
PRINT i,„7 
AOlJ= ; 

A W i .1 “ ^ 

SWL = '.J. - 

IFd.Nc.l )SWL = 1.0-AV(1 )/AV{l) 

AV(I)= 

ILK=:/ 

3ATL=0,0 
SATO=J.D 
DO 3J} KRS=* ,N 

: CALL POINTS IN PROPER SeQUENCE 

IFIKRS.LE.ILK) GO TO :.00 
DO 305 K=KRS,N 
IK=K/JR+l 
JK=K-JR*( IK-1 ) 

IF( JK.LE.0)IK=IK-1 
IF( JK.LE.O) JK=JR 
IF( JK.EU.JRIGU TO 320 

OPE = EP UK t JK )-EP ( IK ,JK+l)+CAP (IK, JK) 
IF(DPE.GE.0.U)30 TO 320 
305 CONTINUE 
320 ILK=K 

MK=ILK-KRS + i. 

DO 350 KR=i,MK 
I J=ILK-KR + ^ 

J = I J 

IFdCP.LT. i.5)GJ TO 2 
PRINT ID,IJ 
2 CONTINUE 
IK.= 1J/JR + 1 
JKsIJ-JR*( IK-1) 

IPiJK.LE.O) IK=!K-1 
IF{ JK.LE,tl)JK=JR 
i ' 1 JA*sI J 3 

VOI« 3 «, 0 . 

BASE«010, 

IF(XCJ,73nEQ.4.0) GO TO 291,,’ 

C COMPUTE FLOW RATES FOR OIL AND WATER 
CALL FLOW '■ ■ ■ 

VOI=X( J,83), -■ E.tE' ' 

YWI =X { J, 84 ) ' ' "iE: " 'E'",-. . r. ,, , - E'-'- 

IFllK.EQ.DVdl^'OiO A; ", '' 


pj P p * f _ A < 

EQUATIONS OPPOluP.i 

-3PP0. jtr.j 
3PPEE 

■ bPPT.' -T 

OPPO' 'i r ' 
JHP0E..9E; 

rjpp 

0 P P u . 

d P 0 , ^ ‘j 

SPPU,'., :> 

SRP^'U. + 

OPP'U. i 

0 PPO> 

3PP0-, ^• 

3PPDlie■ 

LiPPEl, .9 . 

BPP0:;.20C 

3PPfJ.lH.i 

3PPG1.2*: 

BPPiU 2 0. 

BPP0I2AU 

BPPO'^B 

BPPC126 

13PP01270 

BPP01280 

BPPD,»,29. 

BPPai3.> 

DPPOUIO 

3 FP 01 32 0 

BPP01330' 

BPP0i34E 

BPP01350 

3PP01360 

3PP0i37.> 

bPP0138:,> 

BPPOliai 

3PPQ139..i 

BPP01391 

BPP01*-),, 

BPPOIAIO 

3PPLI42S 

BPP01430 

BPP0144J 

BPP0l<^50 

8PP0147j 

8PP0,143D 

3PP01490 

■ ■ BPP0’5'.‘. 

INCOMING OUTGOING 8PP015:. 

dPPOiS/... 

3FP01.3"u 

•3PP.,., 

BPPCEOE 



FORTRAN SOURCE LIST CALC 

SOURC-: STATi f’'. NT 

1F( IK.l, U. j. )VVJ1 = T. 
iF{ iK.LO. ..)u'j r:j .;*** 

UPOAU: K (7A) =LI uoi (jS :..NTERING THE JUNCTION POINT 
I F ( V 0 1 • i. 0 » * .A N U » V WI . £ Q « y , ') ) X ( J j 7 3 ) == 4»*' 

IFIV.jl .=:0. . .AiTD.V^I .NE. J,rj)X( J,73)=2.u 
IFIVDI.NL. . . AND.VWI .NE. U.iOXI J,73)=y,0 
I F( Vt)i .Nt.d.'j.AND.VWI .EQ.O.O)X{ J,73)=1.0 
VOI=A3S(VOI) 

VWI=ABS( VWJ ) 

IFdJ.ST.JlMJGO TO 29.;: 

COflPUT-,,. OUT FLOWING VJLUMcS 
VOO=X( J,9G ) 

VW(J = K( J, ) 

VOO = ABS( VO .i) 

VWO=ABS( VWO) 

DF.LTA=VOI+VWi-^OO-VW0 

IF(V03-LF.. . : -E-IA) VO J=0.3 

iFivwo.LE.n. : 0£-;o) vwy=o.o 

IFIVWI .LE.O. IDE -10) VW 1=0.0 

IFIVFII .LL. . p-Ij' ) VDl = y.u 

IF(ABS(DLLTA).Le. .i _-iO ) DELTA = 0 .0 

DELO=\/UI-VOD 

DELOW= VWI-VWO 

DELO=ABS( DELO) 

DELOW=AB'S( DtLDW) 

IF( DLLO.LT.O.xE-lO) VUD=VOI 
IFfCELOW.LT. .lE-i:) VWO=VWI 
IFi ICP .LT.15)G0 TO i 
838 PRINT37»VOI,VL)0,VWI ,VWO, DELTA 

IFI DELTA. N£. 0. 0 )PRINT37, (X(J,KA),KA=75, 82) 
IFIOELTA.NE.O.OIPRINTO?, (X( J,KA),<A=65,72) 
IFCDELTA.Nc. . .0 )PRINT37, ( X( J,KA+94)tKA=l,32) 
IFIOELTA.EQ.: . ) PRINT37, {X( J,KA) ,KA=75,32) 
IFIDELTA.EQ.O.) PRiNTOT, (X( J,KA),KA=65t72) 

I CDNTINUL 


’ll 


dPPi.. yo 
oppo,. 5*; . 
dPPO jl'' 

B p ° 39 

dP^C' c., J 

OPPO' 61 
3PPL*6.' 
bpp 01 C/-’: 
B PPOi 6-. .. 
3 F P C 6 '+ 

3 P P M D 3 
GPP0.i 3d , 
B POO: i7y 
3PPGi 5:; 
BPPCi d9 
3 PPOI 7 .0 
appoi7io 
3 PPOi. V 2 
BPP'}] 7 3 ^ 
BPPOa 7h-J 
BPPQi 730 
3PPCi76 
BPPO] 77..' 
3PPOi7BO 
3PPOi790 
3PP0.1B..'., 
BPPOMOi 
B PPOI. 810 
3PP0i320 
3PPOi83D 
BPP0184’) 
BPP0l84i 
3PP01342 
3PPL165,,' 


UPDATE OIL CONTENT OF TUBES 
DO 21 K=l,8 
DEN(K)=X( J,K+5S) 

21 X(J,K+56)*V0F( J,K) 

DO 470 JP=I,8 
ICAN(J,JP)=D 
AWDIFI J, UP ) = 0.0 
AOIFOI J, JP)=C.3 

IFU(J,JP + 16).S0.&.0)PRINT lltIJf JP 
',470 CONTINUE 
'‘;C'AS|«0.0 

IFCIK.NE.I.AND. V00.EQ.D.3.AND.VW0-EQ.0.D)X( J,73)=4.a 
I F f * i • AND. VOO* EQ .0.0. AND. VO I • EQ.O .O.AND.VWO.EQ.O 

1EQ.O'4.0')X( Jf73)=4.0 ' ■ • ■ 

IF ('XU'*73).EQ.4.0) CASe=l.O . 

IFIXU'VtSI .eQ.4.0) GO TO 27i..-'; ' •' 

COMPENSATE FOR IMBA ANCE I« INFLOW-OUTFLOW AT A POINT 
MONTE' C'A'iiiLO TECHNIQUE ,■ ' 

I F { I J . Li,,*' J R ) > c A s e = 1 . 0 ■ -- ' '■ 

IF{ IJ.LE.iR.|''GO TO 48'>5' V " 


8PPOl85i 
8PP01852 
3PP01B55 
BPP01357 
BPP0166L) 
BPPOi870 
3PP0188P 
3PP0189^ 
BPPC3 90J 
B PPOI 910 
3PP0192W 
3 PPOI 93.. 

•O.AND.VWI. BPPCl94b 
B PPOI 950 
BPPO! 96 
BPP0i97 

BYPICKING TUBPPOT'^a;’* 
BPP 0,1 990 
3PPC-. 

BP PDF . 



SOURCE STATE f'l MT 


FORTRAN SOURCE LIST CALC 



I Five I .CT .Vji 
IF( Vjl .Z'l.Mj 
IF ( { VOl+VW J ) . 

WDI F=\/wl-yv;j- 
oi Fn=y .jI-v.jl- 
IF( WraF.L!;., . 
iF{ CIFj.Li . 
I F( Cl F J. Lu_ . , 

CALL :iU' 0’:G . 

1 F { D I F 0 . N L . . . 


> A i < U < 
. U '’'i -J > 


VVvI .G- .\/i^O)CAS 
Vl^I.GE.VWO) GO 
LE.0.03 go to as 05 
■( J - LT A*VWI ) / ( VOI+VWI ) 

( rc LTA *Vl;i )/(yoi+ywi ) 

r - , ) WD 1 F = ■ . 0 
5CIFOO.O 


= 1 . 

TO 


u 


^>1 D • W[) I F • L L 


)G0 TO >3. 


w Cl F . NE . 0 . 0 ) P R I NT37 , I X ( J , K A + 7A ) , K 4 = 1 , D ) 

TREaM DATA HAS TO BE CORRECTED DU:: TO FLOW 


1005 
1 00 5 


■ ' • 0 R « 

':.AGS UPS 

805 no 1 ■: 3 K=,,d 

X( J,K+56)=f).. AO 
X( J,K+J 6)=FLU( JfK) 

IF{X{ J ,K + FA) .&£ .1 .5 )GU TO 
IF (X( J,K+6^-) .LF IGD TO 
X ( J,K+tj j )=X( J, K + 94 )+;<{ J,K + 56}-V0F { J,K ) 

X( J,K+i.iC )=^X{ J, K + i02)+yDF{ J,K)-X( J,< + 56) 

X{J,K+56)=V0F{ J,K) 

'L’S CUhTIMUF 
27 V CAT = .. . ' 

CAL I COMDCl ) 

C A T — u • 0 

IFIMMK.FQ., .AND. BASE. LQ.D.O) GO TO 291 
IF(BASc.EQ.l.u) GO TO <t79 
ICAN( J,M)=0 
GO TO 483 

479 00 481 ICA=6,b 

IF( BAS E.EQ. 1.0) PRINT ii 1 1 J » I C AN { J, 1 ) , I CAN { J , 7 ) , I CAN{ J 1 8 ) 

IF(BASE.EQ.1.0)BASE=2.0 

M=ICA 

IF( ICA.EQ.6)M=J, 

MMK=rCAN( J,M) 

ICAN( J,H)=U 

I FtMMK.EQ.t ) GO TO 482 

WDIF=AWDIF ( J,M) . 

OIFO=ADIFO( J»H) 

IF{X( J,M+64) .GT.0.0 )CASE=2.0 


LPP'J,;. , .r 
BPPO. .i. 
d PPG ■ i'9'. 
3 PPG,' 

P P P f. ' . 3 
BPPOt Jh.- 
3 P P 0 ; f A: ■. 
3PP-:,: , 

bPPCt,,-. 

bppo; . 

Cit3PPu> i.' 
b P P G? „ < 

B P P C 2 i 2 J 

3PP0E..dG 

3PP0?’'+. 

3FP02i5 

3PP0I4 6 0 

BPP 0217 G 

BPPO?,B. 

BPPC2.1? 

BPP02l JL 

3 P P 02 2 1 0 

BPP0222 

BPP0223: 

BPP022*tU 

BPP02250 

BPPOZeft 

bPP0227'j 

BPPO22B0 

3PP02290 

3ppo23aa 

BP P 02 310 

BPP023Z0 

3PP0?3'i0 

BPPC23A: 

BPPfiZBS: 

3PP0236U 

3PP0Z36* 


IFf CASE. LT. 2.0) J = MMK 

IFIWDIF.LE.&.O. AND.DIFO.LE.0.0) GO TO 482 

DO 478 K=l,8 

FLOf J»KJ=X(J,K+16) 

VOF( JtiK)«XCJfK+56) 

V iFf Xf U,K-f64) .GT.OIOIAPI J,K) = Xt J,K + 64) 

478 CONTINUE 5 g 

430 IF{WDIF.LE.0.1E-13,) WDIF»O.0 ' 

IF(DIFO.Le.O.lE-l0) DIFO-D'.0^,0 . ' 

IFIDIFO.LE.O.O. AND. WDIF.LE.O.O) GO- TO 485 
CALL CHANGE , ■ ,1/-' 

CENTER CHANGES IN VOL OUTFLOW/ INFLOW, &H0 FLUID CONTENT 
435 CDijTINUE 
48 6 DO 105 K = ! , e 

X{J,K+56)=D£N{K) 

X( J,K+, £.)=FLf,M J,K) 


3PPC237 
BP P 02 38’' 
BPP02340 
3PP024"'.., 
BPP0<f,4l .< 
BPP0242D 
3PP02 43C' 
3PP<)?44..,^ 
BPP0245T 
BPP 02460 
BPP0247.J 
3 PPG? 4o 
BPPC2H-9 
BPPOi. 5CJl! 
oppn^-r. M 
BPP'r’3 - , 
c H p T- ' . : 



FORTRAN SOURCE LIST CmLC 



IF{ X( J,K + C^f) .Gi? . 1 . 5 )GO TO 

IF(X(J|K + 6- ).LE. -.0)Gij TO i05 

/J i J ) =X { J, K + 94 )+;, ( J, K + 56I-V0F { J, K) 

X( S)=X( JtK + „ . 2)+VOF( J,K)-X( J,K+56) 

X { J,K + ‘:j6)=VilF( J,K) 

CDNTlMJi: 


3HFC ' - 
3FPO:o"- 

e p p _ r 7 

3 PPG, 
3PPO. e: 

OPPt;-, o. 
DpPOi cl . 

3PPC; 69 
3 P P C L , 

R poo.' U 

3PP0:;6V 

BPPO/LH 

BPP0;.690 

c»PPO?.70 0 

3 P P ti 2 ( 1. 

BPPC127^ - 

8PP027 y.j 

3PP0Z 7'tO 

3PP0279J 

1BPPCL742 

BPP02743 

3PP0275:.; 

BPPQTTb,' 

BPPC277;- 

BPPO27C0 

dPP0279P 

3PP0280:; 

BPP028.*i 

8PP02810 

3PP02Sii 

3PP0?.8?.:i 

BPPO^Bi-: 

BPP0284U 

8PPC2a5 .. 

BP P €2 6 6^' 

3PP02870 

3PP028eo 

3PPC? 89- 

BPP029:.i '• 

BPP029j 0 

BPP029ii 

3PP029i3 

bPP*1292 

6PP029,;.u 

BPP02 93:. 

B P P 'u. 94 
bPPU>95 
BPP02 9o;/ 


IFtMMK.FQ.O.Ai'jO.BASE.cQ.O.O) GO TO 291 
IF{Bfi,S6.LQ.:.0) 30 TO 479 
ICAN{ J,M)=. 

4B2 CONTnUE 
48 i CONTIMUE 
4S3 BASE=D,0 
MMK=.' 

J=IJ6 ■ 

29J ID=IJ+JR 

DO 330 K=6,9 
IA = K 

IF(K.= y.9) I A = 1 

AV{ I )=AV( I )+X( JtIA+56) 

IF(X(J,IA + 56).LT.0.t.0R.X( J,iA+56J .GT.Xi J,IA+48) )PRINTa.-j JtIA 
1F{X( J, IA+56).LT. O.O.UR.XI J, I A+56) .GT. X ( J , I A+48 ) ) PRl NT37 , X ( J , I A+4 
1) tX£ J, IA + 56) 

330 CONTINUE 

IF{ I.EU.l .AND.ICP.EC. , ) V0LP0R=AV(1) 

C COMPUTE TOTAL OUTFLOW OF FLUIDS FROM THE SAMPLE 

IFIVOI.LE. . i: E-10) V0i=0.0 
IFCVWI .LE.O.IOE-IO) VWi=0.0 
IFdB.LE.N) GO TO 292 
IF(ICP.LT.15)G0 TO 3 
PRINT 37,VOI , VWI 
3 CONTINUE 
AOO=AOO+VUI 
AWO=AWO+VWI 
DO 290 IA=7,9 
K=IA 

IFfK.6«.9I Kal 
SATL*SATL+X( JtK+56) 

SAT0«SAT0+XtJ,X+48) 

290 XOfiTI'NUE 
292 CONTiNi«e''M: 

Cvv;'2 ^PRINT FL0|3 
i##’’ IFflCP.LTil 

'PRINTS?, {X(U,KA I, f<A«57, " 

777 PRINTl2,(XU,KJ,K = i?,32^-V••';:>^; 

IF(X(J,73).EO.4.0I ? PRINT52*I,4:4 I,.. 4 

CDHTINUE " .X,/: 

IF( I.EU.1.4ND.ICP.EQ. DGO TO 2222 y' 

GO TO 79 ^ ’ 

2 222: CONTINUE ■, v-X 


.vr* 




FORTRAN 


SOUR 


IF(X{J,K+;iS).LQ.i,. ) X(J,K + iS) = 25.j 

212 CDUTIMUi:; 

213 IF( K.Ma.r.} GO TO 79 
UG ;r 1A = 7, V 

K=IA 

IF{K.[-0.9 ) K = ., 

lF(X(J,K, + i.6).2Q.I.O) J,K + I6)=3.5 

2lb CONTI NiUt 
79 CONTINUE 
::{j,«s)=o.D 
X ( J , 6 ^«- ) = 0 . 0 
X( j,93 )=i . 

X( J,9£<-)=:% . 

VOI = ,. . i 
VWI=0.J 
19h CUNTlMUt 
350 CONTINJC 
330 CONTI NUf: 

PRINT37, AV( I ) 

SO=AV{I )/VOLPOR 
SW=1.D-S0 
PRINT 50,SO,SW 
69C CONTIMUt 
WI { I)=AWO 
WI { l) = WOa )*ARAT 
CUML=( STQ+STW)/AV( I } 

wo( n = Aoo 

STW=STW+W1 (15 
STO = STO+WO{ I ) 

PRINT 37» W0(1) ,W0{ I ) ,STO,STW,WI { I ) 

w=wo{i 5 *a. . 1 

IF{ I.GT.2.AN0.AQ0.Le.W) GOTO 6TD 
7D0 continue 
800 CONTINUE 
930 CONTINUE 

CALCULATE AND PRINT POROSITY AND PERMEABILITY 

AJR=JR 

AIR=»IR 

PORsAV(ir/{ £ AJR-l. ) =«<( AlR-l. )*GRID**3 5 
P|tW«WIU)’«'VW*C AIR-i. )/£TIHE*GRID*(AJR-l. )>!'CP5 
P||tf.,!.l2r''P0R,PERM 


-3PP03'''.,- 
3 P P ^ ; 
BPPC;' 

cpp--‘;,‘ : 

BPPO--" P-.- 
3 P P C ■' 3 
B P P 0 3 -7 
H P P 'J J :> >. 
BPP03') ^1. 

B P P i* j. i 
BPPv 3„£. 
'dPPo;- . . j 
BPP03, ‘^■L' 
BPPP'-’.S 
BPP03*6C' : 

3PP0317J : 

3PP03itj ' 

bppci:m9. 

BPr-’O^^uO 

3PP032;-0 

BPPQ322 

RPP032^. 

BPP032A-0 

3PP0325^ 

3PP0326, ‘ 

BPP0327: 

BPPO'32uO 

3PP0329.. 

BPP033'3;j 

BPP033i>; 

BPP0332a I 

3PPCB3E. 

BPP033‘+„- 

i3PP03350 

BPP033P0 I 

BPP0337, j 

BPP033a,; I 

BPP03390 I 

3PP0340U I 

BPP03A-1'. ; 

BPPG342:.> S 

12/0 +J 


TBMAP ASSEMBLY CALC 




ERROR 1311 ISN-; 


10 VARIABLE REDEFINED IN ITS OWN RAM 


/ERITY W 


FORTRAN source LIST 


/ / 


SOURCE ST ATtMirSJT 



6 P P ’”1 •+ - 

JR 5 + -j 

rj P O Q _ -l. 

3 p ,j CE ■•;• 7 

3po,/,vy„ 

3 P P ' ; : - -i 
BPP0-. SO 

3 F P .. 3. ,-i 
t! p p :■ .. 
3PPC.'J'- 

3PPvi. A. 

DFP J-'SS 


$ibftcsub : 

SUBROUT IN:- SOLVe 

COMM:t'v|/SYMARC,/X{5.), 16.,J ,AP(50,S} ,GC,VQ,VW,TI>1t,^JF<5 , 
:Xa ?EP{ . j,i0),HYST,CAT,GAMHA,J,FLO(5O,8y,GRID 
CC'HMa4/PKLSS/rjN, A( SO, sO) ,AB(50) 

CaMMO''J/SGL/AA( 5 5' ), 3B ( 5: ) , A V ( 4’' ) tWR ( 40 ) , WOT 4 : } ,WI(4‘ 

COHM3N/SOLB/CAP ( J Sf )» I»iCHANG, ARAT,TRAP 
CQMY30/S0LC/AU3 ,CP, JI ICP 

CUMM3Nl/CHIf|/ICAN( 5a,B ) , A^^ D I F { 50, 8 ) , AD I FD ( 50, 8), BASS 

imtcger count 

.1,. FORMAT (?I 3) 

I ? FORMAT { RFl . 5) 

50 F0RMAT(20X,^5H3IL SATURATI3N=,F10.4,17HWATtR SATURATIC 
75,.' FORMAT (/7'X,*VISC0SI TY RAT I 0( 01 L /WATCR 3 = *, F6. 2 ) 

C FRAME SYST..M OF EQUATIONS 

BASE = 3.-' 

N=IR*JR 
NN=JR=i‘ ( IR-: J 
1S= JR+I 
IT=N-JR 
AV(I )=U.O 
IC=0 

32 CONTINUE 
C0UNT=0 
IC=IC+i 

IFCIC.EQ.II PRINT 763, ARAT 
DO 65 1 JK = :i,N 
IJ=IJK 
J = I J 

C FILLIN PROPER CONDUCTIVITIES 

CAT-IC 

IF( IC.EQ.l 3CALL CONDI,;.) 

DO 65 JJK=1,N 
55 AA( IJK,JJK)=0.0 
195 DO 85 IJ=1,N 
K«I J 

IK=IJ/JR+l , 

JKwIJ-7JR#I IK-l) 

1F^0(^,LE.0> =IK=*IK-1 

;'IFI4K*kf 

C :,R£Aa' 


3 P P 0 7 s 9 

BP P 036 

bPP0.i6„ 

Bppo.'''6':r- 

3POUT6 } 

6PPQ36'+ 

BPP0.i65u 

3PP0366U 

3PPC367 

8PP0368 '3 

BPP03690 

BPP0370U 

3PPQ37 i 

BPP037e... 

8PP03730 

BPP03740 

8PP0375 ■ 

BPP03767 

3PP0B770 

3PP03?Bv 

3PP0379: 

BPP033'.i.' 

8PP033I0 

3PPU3e2.,* 

bPPOBbS • 

BPP0SB4.‘ 

bPP03&50 

B PP QT B6 . 

BP P 03 871/ 

BPPftSBo’ 

BPP03B90 

3PPC3 9-.;, 

b PPU-' 9. 

BPPOr.OTf 

3 PPOE US 

BPPUSV-.' 


CPT«a.i'‘ J,c;: , 

IF( I.EQ.l.AN0.IC.|Q'.i*Ato.|C;^7»ea.l) SO TO 87 
COMPUTE OIL AND WATER SATmAiflON 
I JK=IJ 


( I JK,<2+56) 


AVI I ) = AVI 1 ) + 

io cof-iTiNii':, 

SB CONTINUE 

IF! IC.-O. . ) 
FILL DTLP FO 


„/ / 


SOURCE STATEMtNiT 


FORTRAN SDURCt LIST UB i 


CALL PUIF 

AVOID INFLOW FROM OUTLET END OR FROM INLET END 
UB- S 
K=I J 


DO 76 KL = 2, 5 

I F{ K.GT. JIM) GO TO 76 
CUL= .. 

IJK=IJ 

IF(X(iJK,KL + 2^^).Gt. 1.5) GO TO 76 
As=n.,/ 

IF(K.GT.I';iUAOD. K.LE .J IM) AS = I.r,» 

IF( AS.E0.1.0.AN|D.X( IJK,KL + 64) .LT.0.0) AS = 2. 0 
IF(AS. ^O.L.U.aND.KL .G:-.3.AND.KL.LE.5)CDL = 0.5 
1F(K.L.,.JR.A^/U,X{IJK,KL+6^^).LT-■>.:I ) COL=.;}.5 
IFIXU JK,KL + I6) .LL-.2.U) GO TO 55 
1F{X(IJK,KL+16).GT.4.'.;.)G0 TO 55 
IFIK.Lc. JR) GO TO 55 
UP = 0 . '3 

DP=X(I JKtKL+84) 

OP = DP<' { J..+TRAP) 

ADP = X{ I JK, KL + 64 )-X( IJK,KL + 84) 
IFIDP.LE.U.O.AMD.ADP.LE.O.J) GO TO 55 
IF( AOP.GE.U.O.AND.DP.GE.O.O) GO TO 55 
UP=ABS (DP) 

ADP=ABS(ADP) 

IF( AOP .LT.OP ) COL = : .5 
55 IF(CrjL.F:Q.O. 5) X ( I J K, KL + 24) =1 .5 
IF(COL.FQ.D.p) r.G = 0.5 
76 CONTINUE 

IFJCO.EQ.u.J ) GO TO 77 
CAT=0.0 


COUNT=COUNT+i 
J = I J 

CALL COND( 1 ) 
77 CONTINUE 
B7 CONTINUE 

DO 82 IL=lt8 


IFCXdJ t7i) .eO.4.0 ) GO TO 82 
IFCXil J,lL+24).Ge.l .5 ) GO TO 82 
CPA*X( IJ,lL+64) 

I J f IL + 74) 4CPA 

VS* S.+ Xa; 4 »IL+ 74 ) 

JKJ=CPA 



DO 

X ( rj»i. . ^ . ,,, 

X( IJt<+102)«0.0:7l'’7V: 
X( 1 J,<+< 10 ) = 0.0'^.- 
X{ I J,K + n8 ) = ^.3 
iA=iK '■ , -r 

.JA = JK 

IFfK.CQ.l) iA=IK-l 


■i P P C:> ' 

TRAP APPRUPKIA f. 

■ Hi' P 

appc-u': : 

3 P P J 'f ^ 
ti?P‘ •• ' 
OPPO'vi'i' . 
3PPCi 6 
3 P P V -i- i' 

3 ? P i'-t '■ 

5 P P O'i- D ' > 
SPPt.i''"' ' 

6 PPAi -!■*.!, I 

3PP07,.: 

3PP '.th i A 

fopPOGiS 

3PP0/-16D 

3PPG4..0., 

3PPi,.Ai.7 

3PPG4.P 

BPP04i.9 

3PP043 J 

3PPtJ42l 

BPP0422. 

BPP04230 I 

BPP 0425 . 

BPPOR/Jb ’ 

tiPP0427U 

BPP 04280 

3PP043V - 

BPP 043 l,i 

8PP0452:- 

3PP04330 

3PP0434 ( 

BPP0435; 

aPP 0436 U 

3PP043?- 

BPPG43a ■ 

BPP0439 

8 PP 04400 

3PP044i'., 

BPP 0442 :.,* 

BP»04 43.. 

BPP0444C 

3PPC'445 ' 

BPP‘.44d 

B P P G4*f7 i' 

3PP06 4o.,' 

GPP. .'4 '*9 

L'PP. -7 

GpoO- 

■■'5 ■ 



FORTRAN SOURCE 

LIST U8 * 



SOURCE ST 4 tool NT 






1F(K.L..T) J';=JK+;i 




3Pp'/-:-T- 


IFIK.ST .5 .A.lD.K .L:;, 7) JA = JK-1 




BPP -rSl 


1 F { K , 3 C . 0 . N i_!, K . L t; . 5 ) 1 A = I K + 1 






1F{K.3'.7) lA=rK-i 




rP'-^L'-. 5/ 


I JC = ( I )«J'> + JA 




dPPO-S-/.. 


AA( I J, I JC) =~X n J,K+74 ) 




^ P P 0 > ■ S 


IF{ IJC.L;;,,. JR) V A = \/A +A#. { 1 J , 1 JC )*CP 




7 

CLlNTiNlU : 




BPP.---0. 


BB (I J) =-VA-CPl 




3PP*)~o.: 

35 

4A(IJ,1J)=VS 





3 6 

CUl.lTiNOS 




'3 p P 1' '1 


! F( IC.N" . ■ .AMD. CDUNT. Q.J) GO 

TO 

78 


ri P L 'tt> -j 


00 66 IN = I,'iM 




ap»f>'-6oi. 


INR = l:Ni+JR 




BPPO'-j/v 


AB( IN) =RB( INR) 




3 P P C ^ '3 


DO 66 JN=:l»t^l!i 




b P P L + i-v 9 


JNR=JN+JR 




BP PC' -Yu.:) 


A{ IN, JN)=AA( InR, JNR ) 




PPPOATiU 

b 6 

CONTI NUi 




3PPC4?. 


CALL PRSOL{ X ) 




BPP047:' 


CALL PRSQLC:) 




bPPU-'7-J 


DO 75 K=1,N 




BPP 04751’ 


8B(K)=CP 




3PP04.76J 


IF(K.GT.JIM) bB(K}=U. 




BPPG47?;. 


IFCK.LE.JR.OR.K.GT, JI 0) GO TO 

73 



BPPO^hTSO 


KK=K-JK 




BPP04790 


8B(K)=AB(KK) 




BPPC'+8. .. 

73 

IK=K/JR+l 




BPPCAaiS 


JK=K-JR*( IK-1 ) 




BPP04oi0 


IFUK.tQ.O) IK=IK-1 




BPP0483l 


1F{ JK.EQ.O) JK=JR 




3PP0484G 


EP( IK, JK)=BB{K) 




3PP0485:' 

75 

CONTINUE 




BPP04360 

78 

CONTINUE 




3PP0437 , 


PRINT li), COUNT 




BPP0488P 


IFCCDUNT.NE. J) GO TO 80 




BPP0489> 


IFIIC.EQ.l) GO TO 80 




BPP04900 


RETURN 




3PP0491L 


END 



. 

BPP0492 : 


^ ' IBMAP 

ASSEMBLY 

UB I 

I.Z/U4 

ill^i 

iTX' , / ' : : 






-X Bit ! 





■ES 

P.QR'.'itQyf: -ASSEMBLY 








SOURCE STATCUlMT 


FORTRAN SOURCE LIST 




SIBFTCSUB :> 

SUi-.R iJTi !■!.. Flow 

COMM'-lM/SYnARG/X { S if> .■ ) , AP ( 50 , 8 ) , GC , VO» VW, T I M:. » V OF { 30 
.XA »:.■■) ,HYST, CAT, gamma, J,FLD(5..,B), GRID 

C FLOW CAPICITY OF TUBES 

AHAT = (/LI/VW 


• j = I R* Jk 
IB = N-2=!'JR 
Xi0,03 • 

X ( J , 9 ) = U • 1.' 

DU 2<'j K=: , S. 

VOF{ J,K)=X( J,Ki-56) 

AP( J,< )=". 

IF(X( J,K + 64) .GT.0.0 )AP( J,K)=X( J,K + 64) 

X( J,K+-94)=0.0 
X( J,K+i 2 )= ... 
x{ J,K+iir- )=‘ . 

X{ J,K+.} B )=rj.O 
1F(X( J,?3 ) .K-j.A.C ) GO TO 2 .<0 
C CASF 2 . ' AMS BACK FLOW 

CASE = 3.''' 

FLO{J,K)=X( J,K+i6} 

IF(K.Gc..3.AND.K.Lt . 5)CASE = I.O 
IF (CASE.EQ.1.0. AND. AP { J , K ) . EQ..O.S ) CAS£ = 2.D 
IF (CASE.&Q.Z.vf) API J ,K )=X( J,K + 64) 

1F{X{J,K + 24).G£.?..0)GJ TO 200 
VA=AP{ J,K)4X( J,K+74) 

VA=VA*TIMe 
VA=ABS( VA) 

IFIVA.EQ.O.O )G0 TO 200 
IF(CASE.eu. 2 .», ) GO TO 5 
GO TO 7 

5 IFIXl J,K + 16) .EO.4.2 ) FL0(J,K)=4.2 
IFIXI J,K+16) .EQ.4.1 ) FL0(JrK)=4.1 




IFIXI J,K + 16) .EQ. 3. ,} )FLO( J,K)=3.5 
IFIXI J,K+16) .EQ.3.5 )FLO( J,K)=3.0 
7 IFIFLOI J,K).LE.2.0) GO TO 10 

IFIFLOI J,K).n.4.0.ANU.FL0( J,K}.GE.3.0)GO TO 20 
IFIFLOI J,K),GE. 4.0) GO TO 30 
la IFIFLOIJtKy.EQ.l.O) X{0,K+94)=VA 
IFtFlOI J,K) .EQ.Z.O) X{J,K+102)=VA 

sb' to 

20 IFlFL3U,f<HEQ.3.5) X(J,K + 94)=VA 
IFIFliOI J,t:'K‘yEa.3.3) XI J,K + 102) = VA 
40"tOva.O'.'%As>- 
30,' IFlFt.3,C J#'l^yEQ.4.5) GO TO 40 

IFI,FL0ta*K|.EQ44.l) X(J,K + 94J=VA 
IFIFLOI J#,K74fe0,fc4. a f X( 0,K+102 J»VA ■ 

SO=X( J,K + ^6)/XCi#fK + 48) 

S W ~ 1 . ” So ‘ ' ' '*.('41 ^ -*> , 

AKRW=XtO,K+4:>h'l2-*f$W/ARAT4SW*SW*ll..-2.0/ARAT) ) 
AKRN = XI J,K + 32)*'5a!^SQ':-: ' ' 

AAK = AKRW/VW+AKRH>V 04 rV''C..- • 

A K R W = A K R W 4 X I J , K + 7 4 r 
AKRfl=AKRN4XI J,<+74) /AAK/W':'! 


K p p ! : ■ 

P PPOh. ■ 
),IJ, IR,JK,3PPC:. Y 
SF'P'-vd-'-. 
bPPi-S.!-. 
CPPO 
3 P P C 
RPF'CiPS 
3PPU:'4t,. 

P P 0 -■ + S ^ 

'3 P ^ 

bPPCo-y 

!'.PP05-'s, I 
3 PP Ob 4 9'.'' 
3PP0S3 
bPPCBD... 

B P P 0 3 5 S 

J P P 05 5 i 'S 

OPPCSCY . 

BPPC55:: 

3PD0956..J 

3PP0557U 

3PPC5 53., 

BPP0559 

3PP05600 

3PP056i., 

6PP0362 „• 

BPP0563:i 

BPP05640 

BPPe665..' 

BPPO566'0 

BPP0567D 

3PP05680 

3PP05691. 

BPPC57aD 

BP P 05710 

BPP05720 

BPPC573C> 

BPP0574C 

BPP05750 

3PP0576Q 

3PP05770 

BPPO57B0 

BPP05790 

3PP05S00 

3PPC.58iC 

BPPCsazr.. 

8PP05820 

3PPC583I 

BPPC5 83. 

BPPC5C3S 

BPPOSbOA 

bppobf; 

P.PP'.S ' ■; 


SUHRC:; STATEME^IT 


FORTRAN SOURCE LIST U8 3 




X( J»K+3': )=.' BS(iJiKRvj<5AP( J,K)*TIME ) 

X(JtF+ „'''E)~ABS( AKRW*AP( J,K)*TIME ) 

GD TO 3 

sn=y:iF ( j,K)/K( J,K+Ae) 

sw=l.:'-'';n 

vw ; t=vw 

VI1W = VJ 

swt:T=sw 

SNW=S3 

IF{GC.Lr.u.'J)SWrT = SG 
IF(GC.LT..U ')S'JW=SW 
I F(GC. LT. ) ViMW = VW 
IF(GC. n ) VWFT = VO 
KW;:T = SWrT«S^A:.T ■ 

RNW=... .J^VNW#SNW + SNW#SfiW*{ i.''-Z.D*VNW/VWET} 

TaT= RWET/VWET+RN^^i/VNW 
RAW = R.^tT/VWfc T/TOT 
RAIJ=KMW/VfJW/l UT 
IFISW.L-Q.SWFT) RW = RAW 
IF(SW.EO.SWET) R3=RAN 
lF(SW.fcQ.SNW)RW*RAN 
iF(SW.tQ.SNW) RO=RAW 
X{ JtK + 9A)«R0=^VA 
X ( J,K+i'>2) =RW*VA 
80 VONAX=X { J, K+4a) 

VOMIN^-'l.n 

IF(GC.GT.C.O) VOMAX=w .gO’i'XUfK+A-e) 

IFCGC.LT .0.3 )VOMIN=U. 4*X(JfK+48) 
IFIVDMAX.LT.VOF ( J,K )) VOMAX = X{ Jf K + %8) 
IF(VOMIN.GT.VOFIJ,K }) V0MIN=0.0 
VWMAXsXt J»K+A8) -VOMIN 
VWMIN=X{ J,K+4a) -V3MAX 
IF(FL3{J,K).GT.^.0)GD TO 83 

IF(FL3( J,K),LQ.4.2.GR.FL0( J,K}.EQ.3.S') GO TO 65 
V=VA+VaF { J ,K ) 

IF( V.GT.VOMAX) GO TO B2 
VOF(J,K)=VOF( J, K)+VA 


x( 

X{ J,K+il8)=VA 
IFCGC.GT.0.05FL0( J,K)=4.5 

■'" GO TO 100 
32 CONTINUE. 

XTOfK'FilOlpsV-VOMAX 
; ,''xr'a»k'fii8>?svowAx-voF{ j,K) 

; ',v',#ciFco»::i^T*voMAx 

GO ■TO'ilOO' ■ 

3 5 V=VA+X tj't VOF I Jt K ) 
1F(V.6T.VWHAKI 'GO/ TO., 86 
VOF{J,K^s=VOFt0rR)'-,V.A' 

X( J,K + 11D)™VA. . 

X( j,K+ii8)=o;'0 '' 

GO TO 

36 CONTINUE ■ 

X ( 3 , K-f- 1 } Q ) =VOF ( Jt KT"-VO»lW::','V 





Vs' . 


y'vf 


3PP‘ 

OPPOrr ... 

3PPs c - 

j-i, p p , _ .J 3 

BPP..-: 

3PPGdp7 

dPP!.'..'M ; 

dPPOpJiu.' 

6PP05VJ‘-' 

3pD'''Cy 

OPPtoV . 

BPPOBvA 

B P P 0 ^ V * 

3PPC::9'v 

3PP')59t) 

JPPyiGY 

BPPOBVbp s 

0PPO5 99 

BPP06-, 

8PP06Li) 

BPPOBO^-ij ’ 

BPP06E 

BPP06 

BPP06D5L j 

8PP06070 I 

3PP06f3a:;. 

BPPOB'jV ■ 

BPPOBiOO 

BPP06110 

3PPC6i2v'- 

8PP06130 

BPP06i40 

3PPC6l4.t 

BPP0615„. 

BPP06i6.# 

BPP06i70 

3PPC616.. 

BPP06J.9:. 

BPP062P :* 

3PP0620i 

3PPC6li.3 

3PP0622.,i 

8PPCJ623a 

BPP0624(} 

3PPa6C4j 

BPPa625 J 

BPP 06260 

BPP06270 

BPP0623!.' 

BPPC629-, 

B P P 0 6 3 li F ‘ 
3PP0e5j„L 
3PPy63 ' ' 

B. PP’’::; i 
BPPQ'j v. i 



■o \JI ^ ua M I-* c» -4 0^ Ul - ■ O 0^ \M ii^ rxi * o^ \ J1 ^ U> fv-^ 


FORTRAN SOURCc LIST IJB 

SLiuRCi sTATLMt rr 

X(J,K + ,..:o)=VA-X(J,K+i;-) 

V0F( 

FLO(J,K )=;;..,. 
u 0 T j ^ 

s: v=vor{ j, k )-vgmip^ 

IF(V.Lr.X( J, K + 9A) )S0 TO 84- 
VOF{J»i<.)=VOF( J,K)-X {J,K+-94-) 

GO TO 9* 

34 X( J,K + 94)=V0F( J,K)-VOf IfM 
X( J,K-l- j, 02 ) =V.4-X { J, <+9 't) 

VOF( J, K) = VlIMIN 
FLO( J, K) =2. > 

90 X( JiK+ilO)=X(J,K+94 ) 

X( J,K-Hje)=X( J,K + 102) 

IF{CASt:.?-0.2.. ) GO TO I5v 
;<{ J ,93 )=X( J,9i) +X{ J ,K+94) 

X{ J,9%) = X( J,9*^) +X{ J,K + 102) 

GO TO 200 
15 [i GONTUUc 

X(J,K + 94)=X( J,K + lir. ) 

X( J,K+102)=>:( J,K + II6} 

X(J,83)=X( J,03)+X{ J,K+ilC ) 

X{ J,84 )=X( J,&4) +X( J,K+118) 

IF(IJ.LE.JH) X{ J,K + 24-)=2.;j 
IFdJ.GT.IB) X{ J,K + 24)=2.0 
IF{X{J,K+ 24). £0.2.0) X(J,K+7A)=0.0 
IF(X{ J,K + 24) .tO.2.0 ) X{ J,K+84) = :).3 
IF{X( J,K + 24) .EQ.2.0 ) X{ J , K+24 ) = l .5 

3 X( J,K+110)=0.0 

4 X( J,Kf il8) =!).U 

5 23D CONTINUE 

7 RETURN 

3 END 

I BHAP ASSEMBLY UB 3 


ESSAGES FOR ABOVE ASSE'^BLY 


\OCj 


Si.lURC:. STATcMlNIT 


FORTRAN SQURCl LIST 



$IBFTC$UB i 

SUoRTUTIiJi. CuNO(KKK) 

CUBMD^i/SYH^^RG/K ( 16 ) , AP { 5 G , 8 ) , 3C » VO, VW , FI M t , V3F { 5 ; , 

:XA »FP( ,^-5 ),HYST,CAT,GAMMA,J,FL0(50,8) ,GR1D 
CUMM' iN/CHIn/ 1C\N( 50 ,P ) , AWDI F { 50, 8 ) , AD I FO ( 50, 8 ) , 3 ASc 
GO TO {27 , 79 ) , KKK 

C INTER CHANGES IN OTHER CONNECTED NODE POINTS 

270 CONTINUE 
U FORMAT (POX, El ->) 

13 FORMAT! i/.'Xtcn 3.5) 

N=IR’!'jr< 

JIM=N- JR 
VS=O.D 

IFCRASE.LQ.Z. ,) VS = 1 .2 
DO 26} M==i ,b 

IF(X( J,M+2-U .EQ.l.S ) GO TO 300 
I F( JIM.LT. IJ ) 30 TO 273 
IF(X{J,M + 6-^).LE.9. !) GO TO 28 1 
273 1F{X{ J,M+2't) .GE.2.‘ ) GO TO 230 
300 CONTI NOT 

IF(X{ J,M) .( g.N. i) GO TO 28D 
NX=M + !+ 

IFINX.GT.G ) NX=M-4 
MMK=I J+JR-M+4 
IF( F.EQ.l ) MMK=I J-JR+I 
IF(M.EQ.2 )MMK=I J+1 
IF!M,E0.6 )MMK=I J"i 
IF(M.EQ.7) MMK=IJ-JR-1 
IF(M.EQ.8) MMK=IJ-JR 
IF(MMK.GT.N) GO TO 281 
IFCMMK.LE.O) GO TO 280 
X(MMK, NX+56)=X{ J,M+56) 

X(MMK, NX+1 6) =X( J,M+i6 } 

IF(X( J,M+16) .EO.4.1 ) X! MMK,NX+16)=4.1 
IF(X( J,M+16) .tO.4.2 ) X{MMK,NX+16)=4.2 
IF(X(J,M+I6),E0.3.5) X(MMK,NX+16) =3.0 
IF(X{ J,M+16) .EQ.3.0 ) X{MMK,NX+16J = 3.5 
IF(X( J,73 ) .tQ.4.!)) X( MMK,NX + 24)=2.0 
IF{X{ J,73).t:Q.A-.0) X{ J,M + 24)=2.0 
IF(X(MMK,7'J) .tQ.4.D) GO TO 275 
iFCX{J,M+24) .eQ.1.5 ) X(MMK,NX+24)=2.0 
!FIX( J,M+24) .EQ.1.5 ) X(J,M+24) =2.0 
IFa(MMK,NX+24) .EQ.2.0) X ( MMK, NX+74) =0.0 

X(MMK,NX+84)=0.0 

’ >‘vA:"lF{X(M.MK,NX+24) .EQ.2.a) X( J , M+74) =0. D 
' X( J,M+84) = 0.0 

IFfX(^MK,NX+24) ,EQ.2.0) GD TO 275 
27 4 IFte'ft'T.Ne.l.O) GO'TO 275 

IF(XU»H424) .£0.1.5 )X{MMK,NX+24I=2-0 
I F ( X ( J V»V'2'4) . EQ . 1 . 5 IX < J.f M+24J =2*0 
IF(VS.Ne.O.O) 30 TO ZIS’'' . 

IF{M.EQ.l.0!l*«.GE.7) GO TtY2t2 
X(MMK,83) = X(MMK,83) +X{4,*M4 t|.10,1 
X { MMK, 84 ) *X( MMK , 84 ) ’ : 

GO TO 275 \ y-'' '—'y ■ 


0 P P 0 A "j 
3PP';F5 7 
) , I J , ! R , Jli , b P P 0" i S 

PPP0639:.. ! 

3PP':-jo ; 

3 p P •- 6 6 . 
BPPGSj„ 
HpPOSfcE.; 

3 P P C o G ■■■■ 

3PPG665 ■ 

3PP0666-.> 

OPPOG.'V,/ 

3PP066i 

BPP'OEdY 

aPPO670b . 

BPPONV.t' 

aPPC673 

BPP0674. ; 

dPP067 5.- 

3 PP 0 : 7 S ■ 

3PPA677 

BPPOGTe 

3PP0679 

3PP068.,. 

3PP0681 

BPP0682 

BPP06abc 

3PP0584 

3PPt}6a5. 

BPP0686. 

8PP06870 

3PP068S 

3PP06899 

BPP069:'. 

BPP06910 

3PP0692 ■ 

BPP0693 

BPP0694-;j 

BPP0695'- 

BPP06960 

BPPG5970 

BPPQ69S0 

3PP06990 

3PPO7y0j 

BPP07'Ji0 

6PP07020 

3PP0703C 

BPPC704.‘ 

8PPC7059 

BPP07060 

BPPC7-. 7. 

BPP07 .<3 ' 

BPP07:9 ’ 

3 p p 0 : , . ' 

i H r - ■ , 


no 


SOURCE ? 

>TAT5 

272 CuUTJNU- 


I F { X ( J , F 

' + 1 CL 5 

1 F ( X ( J , . 

4^. 

IF{X{ J,F 

4+::; p 


FORTRAN SOURCE LIST US 


:?/ -/ 


. 10 . 3 . .UK.X( J, M+r 6 ) . EQ.J .5 J GO TO ^'75 
) .NI.K( RMK,NX + 94 n AD I FD { J , M } =K ( J , M+i ] )-X( 


!4MK , */.+ :>' 
MX + 


IF(ADIFO{ J,M).Nr.L.U.AND,AwDIF(J,M).Nfc.C.S) BASi: = i.J 

X(MMK,N/, + 9 h-)=X{ ) 

X(MPKfNX+„ S ) =X ( J, M+i , 8 ) 

IF( BAS F 0 , 1 ' .u) iC<\tN:( J,H) =MMK 

I F ( 3 A S i.'. • 1 . U • ~ • ) F R 1 0 T; » X { J » M+ 74 ) » K { J i M + 64 5 > X { J » M + 1 6 ) j X ( riMK j NX + 7 


.,X (MMK, :s|y + S-> ) ,X( HMK, NX + lb ) 

I F { B A S E . t Q . 0 , 0 . A N D . X ( J » M + 2 4 ) 


I F ( BAS J: 


i.OIPRINT 


INK, ICA 


) ,X(<J,M + 94 ) jXIJjM + iuE ) 


0 . AND.XC J,M+ 24 ) .GE.Z.OJ PRINT i 3 ,X( J,M+i. 


N ( J, MJ 
) , X { J , M+ ,, 


CONTINUE 
CONTINUE 
GO TO 500 

FILL PROPr.K CLN DUG T I V i T I ES 5 

OiJ 78 K = i , <j 

BASE=0.0 

IF{X{ J,K) y.0.0) GO TO 78 
IFtCAT.FQ.j. )X( J,K+24)=i.3 
IF{X( J ,73 ) .t.;u.4..') X( J,73)=3.0 
SO=X{ J,K+56)/X( J,K+48) 

SW=1.0-S0 

IF(X( J,K+16) .tQ.1.0 ) X( J,K+74)=X( J,K+32) 
IF(X( J,K + I6) .EO.Z.'} ) X( J,K + 74}=X( J,K + 40) 
IF(X( J,K+i6) .EQ.I.O.OR.XI J,K+16).EQ.2.0) 33 
IF{X(J,K+16).NE.4.6) GO TO 821 




IF(X( J,K + 16 ) .EO.Z.'T ) X( J,K + 74 }=X{ J,K+ 40 ) 
IF(X(J,K+i 6 ) .EQ.I.O.UK.XI J,K+ 16 ).EQ. 2 . 0 ) 33 TO 830 
IF{X(J,K+ 16 ).NE. 4 . 6 ) GO TO 821 
VL = VQ 

IF< VW.GT.VL) VL=VW 

X( J,K+ 74 ) = X( J, <+32 ) *\/ 0 /( VD- 4 'S 0 +VW 4 SW+ 6 . 0 *VL*X{ J ,K} ) 
GO TO 830 

IF(XtJ,K + I 6 ) .E 0 . 3 .: . 0 R.X{J,K + 16 ).EQ. 3 . 5 ) GO TO 822 
GO TO 823 

X{ J,K+ 74 )=X{ J,<+ 32 ) *V 0 /( V 34 S 0 +VWf'SW) 

GO TO 830 
ARAT=VO/VW 

IFCGC.GE.O.O ) GO TO 824 
IF(GC.LT.O.O) GO TO 825 
AKRM«SW»SW*X( J, K+ 43 ) 

AKRNaXI J,K + 32)4 (2.0>5<SU’«'ARAT+S04S0*{1.0-2.a’i'ARAT) I 

TOT=AKRW/VW+AKRN/VO 

X(^NK+ 74 ) = X{ J,<+ 32 ) WO*TOT 

6 O'TO -830 

AKRW^XC j,K+ 4 :H 4 ( 2 . 0 *SW/ARAT+S 04 SO*( 1 . 0 - 2 . 0 /ARAT} ) 
AKRN^XU»K + 32 ) 4 SO«SO 
TOT«AKBW'/VW+AKRM/VO . ' 

X( J,K+f 4 Ts*t 0 T=«'VO*X< 0 ,K+ 32 ) , ■■'■ 7 , 

CONTIMUe ’ " , 

X( J,K+ 74 )=ABS(K( ■ 




X( J,K + ilO)=O.U , :■ ' ' 

X { J , K + 1 1 8 ) = 0 . a .. :iT,f . A A,-. . z,. . 

X ( J , K 1 84 ) . a T :■ 

X{J,S 3 )= 0 » 0 ' A.- :.T, ■ ■ 


OPPC?’:. Tv 
3PP(.7'' : 

riPP^'? 4 

C P H i . T' !, 7 
3 PP 07 .., ... 
3 1 '.'r’ nT ' 
,BPP T 7 , 
tiPPtY^. 
SPPO '/ 2 .. 
JiKPOl / • 
3 PP *.,7 -+ ■ 
rtpPC 7/_5 
sppo?;;--.,; 
3 PPa 7,-7 
rtFP(j 7 :i .8 
BP POT 2 9 
dPPOTSu'' 
6 PPP 7 ..; I 
E PPG ? 37 
BPP 07 :vu:. 
BPPOTS^s-O 
BPP 0735 .. 

■ BP P 07 it),'; 
6 pp 07 770 
BPP 0738 C 
BPPC 739 • 
BPP074:,iT 
BPP 07410 
3 PP 074 ^. 
BPPC 743 ' 
BPPC 744 ,> 
BPP 07450 
3 PP 0746 . 
3 PPQ 747 J 
BPP 074 S: 
BPP 07490 
3 PP 075 .,}(. 
8 PP 075 ii‘ 
3 PP 0752 'T 
BPP 0753 U 
3 PP 07510 . 
8 PP 0754'.3 
BPP 07550 
6 PP 07560 
3 PP 0757 ':. 
BPP 0757 i 
BPP 07 53 3 
3 PP 07590 
BPPCT 6 -1 
RPPC 76 ' 

■ BPPCTo.. 
3 FP CT 




ill 


FORTRAN SOURCE LIST U6 % 

*7/ -t 


SOURCE STATt ME'^T 



X( J,&4)=0. j 

3PP0:“,- 


iK=j/jR+: 

RHP 07 6J 

;3 

JK=J-JR*{IK-.,) 


•A 

IF( JK.LE.i nK=IK-„ 

3PP076r.' : 

'7 

IFIJK.Ll'.O) JK = JR 

BPF.'7 - . , 

■f'k 

d. 

IF { I K. ■ 0. ) b(i Tu . 41. 

B P P . 7 6 a 

5 

1F( IK..-0.2..J;i:!.K. ;u.l )Gn TO l40 

3P®l! 7 7 i 


IF ( Il< . K.G!- .7 )&□ TO 140 

3PPO/7 

■3 

R=oc*( { .j/:; { j,K) ) ) 

BPPCT^- 

■4 

R = -R 

3 p V r/ 3 

■5 

b*j T:.,J * 

3ppo77‘.-v; : 

■6 (40 

*^' = GC/>< { Jt K) 

Bppczrac ’ 

•7 „3C 

CUNT3 NO:- 

3 P P L 7 '/ L 

:•» 

1F(X{ JtK+j. 6) .LQ.O. } ,;{JjK + 84}= R 

ypDr777 

3 

IF(XU,K+^6).tQ.3. 5) X{J,K + 84) = -R 

BPP077 i. 

6 7a 

CONTINUE 

3PPn7V9u 

D 500 

RtTURN 

3PP<;7h 

1 

lND 

3PPC7i>i : 


I3MAP ASSEMBLY UB 4 


ESS4G5S 

FOR ABOVE ASSi-MBLY 

’ '■ i 

: 1 


FORTRAN SOURC 


LIST 


M 


sour: 


STATi^Mc'IT 


SIBFTCSUB 5 

SUBROUTINE CHAM&E 

COMHDN/SYMfiRG/X(50t 160) ,AP(50,8), 

IXA ,FP{ ;j,i'))tHYST,CAT,GAMMA, J,FLD{5:;,8),o 


59 


6 80 


C 

c 


55 


u 

3 

6 

7 

0 

1 

L 

7 

2 

5 

6 
7 
1 
2 

5 

6 
7 
2 
A 

5 

6 

7 

2 

3 

4 
|5 


55 5 
556 


COMMO^/CHA?,/r.IF :.,WD IF 

' - "• • iO),I, ICHANG, ARAT,TRAP 

1 f ». ft T r* n 


3C,V0,VW,TIMt, YOFCSOtb), 
RID 


I J , I R , J 


iICP 


COMMON/ S0L8/:a? ( . U 
CUi'1 M 0 J / SO L . / i I JO , C P , J 1 M 
OlOEMSUO 0(9),C(fc) 

R.AL K 

iSHuI F0=, El; . 5 , 5X , 5HWD I F = , E 13 • 5 ) 
AND.I.GE. .7. AND.J.EQ.22)PRINT68 
AND. I .GL...7.AND. J.EQ.22)PRINT63 
NO. J. £0.22) PR I NTS a 


A 


FO>,M.;.T ( 

lF-( i:.p. 

I F( ICP. 

1 F (I : p . 

IF( ICP.FQ 
FORMAT { iX, 


^ 0 , 
■U . 
= Q. 


. A N 


t 01 FO » WD I F 
<t (X< J,K) ,K=i75 
(X( J,K) ,K=75, 


DO 


5j'? 


.:.Ai,D.I.GE.. 7, _ . , _ __ . 

?. ANO . I . G E. 17 . A ND. J .FQ . 2 2) PRr>lT6S0 , { X { J , K + 94 ) , K= 1 
FL:*5. 7) 

,3 

FOR FLOW CHANGES SO 
IN EACH NODE OF THE 


55 4 


555 


K = 

SELECT TUBES 
IS SATISFIr.D 
C(K)=D,0 
CON liNUE 

ETTIN3 FLUID IMBIBITION 
BASE=0.0 
KA=C.3; 

CAT=0. 

AC E = 0 • >1 

IF(GC.GE.O.O.AND.WDIF.GT.O.G) 
IFIGC.LT.' ,,*.ANLi.D! FO.GT.D.O) 
IF(KA.tg..t.O )KA = 4.5 
IF(KA,EQ.4.5) 3AT=0.3 
B( l)=D..O 
DO 555 K=1,S 
V=L.D 

IF(X{J,K + 24) .GS 
IF(C{X).NL.C).0) 

IF(FLO{ J,K).Nb 
IF(AP( J,K) .Li:. 

V = AP{ J,K)*; { J,<+7 
V=ABS( V ) 

B(K+',) = r{ K 3 + V 
COOT! aj- 
1F(B(9 ) .EQ.O.O) 

rn=rndy:.. {>0 . ) 

DO 557 K=^,9 

, '.IF{RN.L^.B(K)/3(9) ) GOTO 558 

$5T;C;;0NTINUE 
55 8''K'eiKH 

CrKOt.i.0 
GO T0:,57O : 

PENOUL^R'TLOW CHANGE 
IF(CAT.Gff.D.33 GO TO 560 
KA=4.5 
CAT = .'.3 
GO TO 554 

IFtCAT.EQ.O.S) GO TO 
IFICAT.EQ. ■'.4} GO TO 


THAT MATERIAL 
MODEL 


BALANCE FOR EACH 


KA=1 .0 
KA=2.D 


.2.0 )Gu 
GO TO 
GO 
GO 


. KA ) 


TO 

355 

TO 555 
To 555 


) 


GO TO 559 


5 5 9 


5S 


.561 

562 


Lipp rr. ■ 
3PPC7 ; , 
,,3PPC1 ■+ - 
JPPo7.;i.- 
B P P f 1 . . ;; 
3PP0V -iS; 
li P P 'i i 5 i 
3 P P « < J Y 
ft PPG? > . 

3 PPG? 1 9 J 
3PPr., /3 9. 
BPPC7 j9^ 
BPPr:7o9. 
3PP07 6 9'+ 
3 PP 07 89 3 
BPPO/y 
PHBPP0f9.u 
8PP079EC 
8PPCJ79-:. 
BPPG7V4 
BPP0795U 
BPP079OU 
3 HP 0797 
BPPG795 
BPP07 9B1 
3PP0799U 
3PPG8 j. , 
BPPCSOi ,■ 
3PPOfiO2 0 
3PPCb.^4 . 
8PPC8D4,> 
BPP08:>5 ■' 
3PP0S06G 
3PP08., 7 
BPP08..8 . 
B P P IJ8 ; 8 i 
BPP08090 
3PPC8i 
BPP06i.L' 
3PPC8127 
3PP0Ba3U 
BPPfidlA ' 
BPP€8i5'' 
BPP0316U 
BPP08:.70 
3PPCJB.< .30 
BPPCi319'w 
8PP08200 
3PP062:i0 
3 PPC-8 
BP p Oft/:? 
SPP0bk<;''-0 
SPPLo' u 
Bpi-'-, < ■ 



113 


SQURCi: STATtitMT 


FORTRAN SOURCE LIST UB 5 



IFJCAT.tQ. .5) GO TO 563 
IF(CAr.LQ.C.6) GO TO 564 
IF(CAT.FQ.u.B)3U TO 565 
GO TO b?J} 

561 COMTIMUE 

C SLUG FLOW CHA'-IG,; 

IF(WOIF.GT..:i.,OKi!=4.2 
IFCDIFQ.GT. ) KA = 4.:L 
LAT= 

T.J BS-.' 

5 Si. iF(W01F.GT.O.O)KA = 4.l 
IFCDIF.j.GT. i..JKA = 4.i 
CAT= .. J 
GU TO 

562 CO^iT^■^lU2 

c displaclmlnt flow change 

IFIDIFU.GT.O. j) KA=3. 

IF{ WOI F.GT.G.O) KA=3.3 
CAT = :J,5 
GO TO 554 

563 CONTINUE 

C NON WETTING FLUID IMBIBITION 

KA=0.D 

IFCCIFU.GT.Q ..J)KA = 2 .. 

IF(WDIF.GT..J.G)KA = i .0 
CAT = :i.6 
GO TO 554 

564 CONTINUE 

IF(WDI F.GT.O.O)KA=3 .0 
IF{ DIFO.GT.fi. u)KA = 3. 5 
CAT=a.7 

IF{KA.fcQ,0.0)GO TO 620 
GO TO 554 

570 V=AP( J,K)=«=>U J,< + 7^) *TIM£ 

FLO( J, K)=KA 

IF(KA.2Q.3 .OIFLOI J, K)=3.0 
IFIKA. ,0.3.5 )FLQ{ J,K)-3.5 
IF{KA. :0.4..' )FLQ( J, K)=4.1 
iF{KA.2 0.«t.2)FL0{ J, K5-4.2 
IFCKA. FLCJ{ J ,K )=i.O 

IFJKA. _Q.2.u) FLQ( J tK. )=3.5 

IFtCAT.iiU. . 5 )FLO ( J ,K) = 4.5 
VGMAX = X{ J,K+‘ic5} 

.2. VOMIN=U.O 

2;i23''''IF{GC. LT. \/OMIN=^ .04=!‘VDMAX 

't' ' i'FieC.GE.''. 5 VQMAX=2.96*V0MAX 

iFiXCUf K+56) .GT.VOMAX} VOMAX=X{ J iK+48 J 
lFCXC4#K+56) .LT .\/OMIN) VaMIN=0,.0 
VWMAXi=XI J»K+48)-V0MIN 
VWMIN*X( J»K+48) -VOMAX 
VQLO = VOMAX-X ( J, K + 56 ) 

VOLW=X ( J, K+56)-V0MI N 
CASE = 3.0 

IFiDIFQ.GT. .}.:)) V=X(J,K+102I 
IF(WDIF.GT.Ej. 0) V=X{JfK+94) 


3 p D rjc ■' ; 

0 -^PO'- ; . 

3 P P ' 

BPP'O'-'O 
t'.PPiO', ' 

& PPOf i 
BPPL:'.'. : 

BPPiF 
BPPOc . 4 '• ' 
BpoCt-^r 1 
3PPo:“5:„ I 

B P P 0 . 6 , 

3 

UPPO i.-j5 
BPPueOT,- i 
6PP0.,!J0 ! 
3PP0R390 [ 
3PPCb4’. I 
BPP0e4* I 
3PP0S<+2u I 
3PP0B4- i 
3PPfi&4 + .. ! 
3PPCt45 
BPP0846^i 
BPP(J847 . 
bPPC849 . 
BPP0849i I 
3PPOP‘tO^- i 
3PP08493 
8PP08494 
3PP08495 
BPPC8496 
BPPC85-.,, 
BPPfJBSI ^ ; 
&PP0B52U i 
3PP0853 V 
BPP0854/ 
BPP03550 
BPPO05BO ‘ 
3PPC357<; 
&PP0S5 7.. 
BPPOaSUO 
BPP06590 
3 PP 08500 
BPPCiB6i.7 
6PP08620 
BPP 08630 
3PP0664G 
BPP 08 65 2 
3 P P 0 o 6 6 0 
3PP0267G ' 
BPP086^. i 

Bpr'. J- I 

juoe' ; j" : 



1/4 


SOURCE STATEME’n 


FORTRAN SOURCE LIST UB 3 



IFIWDIF.GT.C.O) G3 TO 601 
IF( DIFU.GT.0.0) GO TU 603 
C WATER CXClSS 

601 IFtV.GT.WDlF ) GO TO 6 2 
WDIF=WUIF-V 

x( Jt <+9 +)= .. : 

X(J,K-(-: ': ) -' { J, K+i ; 2) +V 

IFtV.LT.VOLW) GO TO I, 

FL0(J,K) = 2....‘ 

VOF ( J, K,) = VuH’N 
Glj TO 

1.. VliF( J, K )=X( J ,K+ 36)- V 
GO TD .. E 

60 C VOF { J, K)=VUMIN 

IF(VOLW.GT.WDIF) VO F { J, K ) =X ( J f K+ 56 ) -WDI F 
IFIVOLW.LT.WOIF) FLO(J,K)= 2.0 

IFtGC.LT.O.D.ANO.FLGf J,K) .EQ.2.0) FLO( J,K) = ^+.5 
X( J,K+94) =X{ J, <+94) -WuIF 
X( J,K + 1 ■2 )=X{JtK + 1)2)+WDIF 
WDIF = 3.') 

12 CONTI 'lUe 

SO=VOF( J,K)/X{ J,K+48) 

IF(FL0( J,K).fU. 4.5.AND.GC.GT.0.'>) CASE=C.5 
IFfCASE.EQ.o.S. ANO.SO.LT.0.80 )FLO{ J,K )=4.1 
1F{ FLO{ J,K}.c.Q.-+.2. CR . FLO { J, K ) . EQ. 4. 1 > BASE = 1.0 
IF(BASe.fcU.l..,. AND.GC.GE.t .0) CASE=1.D 
IF(BASE.EQ.l.w.AND.GC.LT.u.C') CASE=2.0 
1F(CASE.EQ.*.0. AND.SO.GE.0.8) FL0(J,K>=4.5 
IF{CASE.hQ.2.0.AND.SO.LE.0.2) FL0(J,K)=4.5 
CASE=}.0 
BASe=3.:> 

IF{ WDIF.NE.0.0) GO TO p54 
GO TO 590 
C OIL EXCESS 

603 IF{ V.GT.DIFO) GO TO 6 4 
DIFO=DIFO-V 
XT J,K + 94) = :;( J, <+9-:f) +v 
x{ j,K+i, r, ) = , 

IFIV.LT.VOL'O GO TO I.: 

VDF ( J, K) = VOf-U-,X 
FLO(J,K}=:. . 

IFCGC.GT. . ) FLn( J ,K)=4.3 
GO TO 554 

13 yOF{J,K)=X(J,K+56)+V 

' '-M ' Gn TO ' i-- 

634 VC3FU»K)=Vl>MAX 

IFlfpLO^GT.DlFO ) VOF( J , K ) =X ( J , K+ 56 )+OI FO 
IFtVOLO.LT.DlFO ) FLC{JtK)=l ,0 

I F ( GG- GE • t> .t . AN D. FL 0( J, K ) .EQ. 1. U ) FLO { J, K ) =4. 5 
X{ J,K+102)=X( J,K + 1.02)-DIFO 
X{ J,K+94)=X( J,K+94) +DiFO 
DIFO=0.0 
1,4 CONTINU: 

SO=VOF( J,K)/X(J,K+48J 

1 F( FL 3 ( J,K ).EQ. 4 . 5 . AND.GC.LE.O. 0 )CASE= 0.5 


3 P ; J 

OFPvll.-- , 

3 p O \ 

rOPPOl ( 3 . 

3 P P C ? i : 
rt P P 71 
ri P P ’3 r: '/ . i 
ii P P 0 . I 9 ‘J 
3 PP l; . ; 

BPPi.'.'i,. ' 

3 P P Oo , 0 ; 

PPPOPtAfj ; 
OPPGGoC ‘ 
h P P B o ‘ ■' 
GPP 0Gb ?o i 
OPPGblo'i ' 
B P P i* 4 C 9 ' I 
BPPOBP - ' ’ 
BPPOBP.,^. : 
3 PP 06 9X 0 !; 
3PP0893. ' 

riPP0R94 ' 

'6 PP 06950 

BPP0&96 .i' 

8PP0897;. 

bPP5398 i 

8PP08990 ; 

3FP€9:y 

3PP(J9bi..., 

BPP0902C 

3PP09030 

3PPU9;>4 

3PPU9 .o..' 

BPP0906U 

BPP09‘)70 

BPPG9;-di' 

bPP09 9- 

BPPOUiOC 

BPP091 lO 

BPPi39I2 ^ 

bppy9j,,3 0 

8PP09J 40 

3 PP 091 30 

3PPU9i.6ri 

BPP09i7>..' 

BPP091B0 

3PPG9i9o 

BPPC92‘i-,i 

BPP09 2j.,.,< 

3PP092^ u 

f)PPC9.>3 

hPP''39>-* 

o P P ‘ ‘ V 

;j .3 fOO ’ , . ■ 



SOURCE STATEMc^T 


FORTRAN SOURCE LIST Ub 5 


’■ 


IF {CASE. no. 0.5. 4N0. SO.GT. 3,2) FLO ( J »K)=4.2 
I F { FL 3 ( J , K ) . f Q . 4 . 2 . OR . FLO { J, K ) . £Q. 4. U 6 ASE = - . : 
IF(BASF.£-U._.i . AMD. GC.LT.O.O) CASE=1.* 
IFIBASe.FQ.l. . AND. GC .GE.3.l>)CASE = 2.3 
IFCCASr.EQ.,^ .0. a,nd. SD,LE.0.2) FL0(J,K)=4.5 
1F{Ch£L.::o.;-....*. 4U, SO.GE.0.8) FL0{J,K)=4.5 
BASc= . 

C fi r* 'r* __ 

H Z - . * 

1F( DiF j. >!..•.). J ) GO TO 554 
CLiNTIMU.,., 

6'^ IFIDIFM.Ll. )D.Fa= 

iFlWDIF.Li-:., .i 't-i )wi;aF=;. 
lF{CIF3.::0.D.0.AM0.WUiF,i-U.0.0)G0 TO 633 
A C c — A C + 1. . 

DO 35 K=i , B 
35 :{K)=D. 

GAT™ 0.3 
K A=4. 5 

IFSACE.LQ.i. )G0 TO 554 

63 I F { WOl F . NE- . 3 , 0 . OR . D IF J, NE . 0. 0 ) PR I .MT690 , DI FO , rJD I F 
RETURM 
END 

FORTRAN PROGRAM UB 5 

i FUNCTION SUBPROGRAM REFERENCES 

f 

i ■ . . 

! IBMAP ASSEMBLY UB 5 


sppLO.. : 

3PP ^9. 3 J 1 
3Pi’09.- „ ' 

3PP . ' 

B P P 9 0 
BPP09;--f.. 
appoo-i^i' 
3ppi;F5B 
BPP09..7.. ’ 
PpP09.;...o 
Bppooio; ; 
3PP4 535:; ' 
BPPCQjB;. ■ 
SPP09>Jh 
3 P P fN J P 5 ! 
3PPC9386 i' 
BPpg933? [ 
BPP 09590 ! 
3PP094.-. ( 
BPPG94:' 

12/.' 

, ' , ■ , ' ' 'I 


1 i/. - 


ieSSAGES FOR ABOVE ASSEMBLY 



SOURCE STATtMnMT 


FORTRAN SOURCE LIST 


, / - 


.1 3FTCSUB 6 

SUBRUUTLMb PDIF 

COMMC^/SYMARG/X (S./, i6 ) , AP ( 53 , 8 ) , GC, VO, YW , T 1 HE , VUF { f 
,'■:?{ j.0,n) ,hYST,CAT,GAMMA, J, FL0(50,S) ,G^!D 
PPcSS UiFFcRCNT lALS FOR ALL TUBES AT EACH ''iOOt PUIOT 


f ci ) , I J 


iXA 
FILL 
J = IJ 
IK= I J/ JR + 
JK= IJ-J< 


i F ( JK. L 


(IK-: ! 

) IK = I K-. 

IFCJK.L-. )J<=JR 

f:jph;t { X, ..1 ) 

>.( J,7A)=i.P {IK, JK) 

1 F { I K. . u . ;■ ) GU To 3 7 
IF( JK. .U.l )G!i TO Pi, 

1F( IK. ..Q.IR} 03 TO 96 
IF( JK. -,0. JK) GO TO iO ,.' 

X{ J,65)=FP(1K, JK}-:-P( IK 
X { J , 66 ) = t: P ( I K , J K ) -F P ( i K, j 
X{ J,67 j=fiP(lK, J!<)-t;P( iK + i, 
X( J,63)=EP(IK, JK )--; P{IK + 1, 


A, JK+ 1 ) 
JK + i) 
JK +1 ) 


3 PP'-V'f.. 

hPP-.9‘'-. 

I K , J 7 , I-J P P C9 -r J J 

^ P P 0 9F 
3 P P ■ I 

b pPfiP^iV ' 

^PPfP' ..J 
3 pp rjOA ;i., 
3 PP-'K-»:j 
HPPiOX. 

3 P P O’-j 6 

3 b p ij 9 ^ 

OPPbP’sp-'i 
:i P P J 3 3 C ■ 
?PPC756j 
OPP 0937 

OPP 0 O.S 9 
3 PP096a>l 
JPPUYo. 


X( J,69)=EP(IK, JKI-L'PC IK+i, JK-. ) 
X{ J,73 )=eP(IK, JK)-FP( IK, JK-1 ) 

X( J,7i )=EP(IK, JK)-E P(IK-i, JK-1 ) 
X{J,72)=EP(iK, JK)-EP{IK-*, JK) 

GO TO l :2 

37 IF( JK.EQ.l ) GU TO 88 
TFiJK.EQ.JR) uO Tu 9n 
X(J,65)=0.(» 

X( J,66)=_ . 

X(J,70)=O. 

X(J,7l)=0.0 

X(J,72)=iJ«* 

X( J,67)=EP(IK, JK)-EP( iK+l,JK+i ) 
X(J,63)=LP(IK, JK)-i: P(iK + :., JK) 

X( J,69} = CP( IK, JK)--:P( iK + l, JK-i. ) 
GO TO i. " L 

83 X{ J,57) = .P( iK,JK)-";p( iK + l, JK + 1) 
DU 39 K=7, . 
iFIK.LO,.:) GO TO 39 
X ( J , K + o A ) = . 

39 C(3fJTi']lJ'i 
GO TO >2 

90 X{ J,b9)=uP(IK, JK)-FP( IK+i, JK-1) 

. DO 9i K = }. , n 

a/V'JF{K.EQ. 5 ) Cn TU 91 
X{J,K+69)— . 

91 CONTINUE 
GO TO. 132 

9:. ifcik.eq.ir) go to 99 

X( J,65)=gP(IK, JK)-EP(1K-1, JK+l) 
X( J ,65 )=EP( IK, JK)-EP(IK, JK+1*. 

X ( J , 67 ) =E P (I K, J K ) -E P( IK+i , 

00 93 K=^^, 3 , . „ 

X ( J , K+6‘'-!. ) =0. 0 '' ; 

9.3 continue . ■ ' . 


3 P 0X96 A • i 

aPPC96b ; 

6PP0y6+U 

3PPC965 * 

8PP<J966 > 

8PP09670 

BPP09680 

3PP09690 

BPP097 f > 

BPP097i0 

BPP09720 

BPP0973U 

BPP097 

3PPa97b-, 

3PP0976: 

BPP0977} 

8PP€978:J 

BPP09 790 

BPP09800 

BPPC981y 

BPPCJ9823 

BPP09830 

3PP098't0 

8 PP 098 50 

BPP0986i 

BPP09870 

3PP09SBW 

BPP09890 

Bppcyg’,.' 

3PP099i0 

3PPC992.. 

8PP6S993', 

6PP099v- 

BPPOOGl " 

J Ppi 



1/7 




FORTRAN 

SOURCE LIST UH 6 

- 1/ 

>N 


SOURCE STATEML'IT 





GO TO *32 


3 P P j ■ 


94 

X{ J,65)=fPaK, JK)-rF( iK-lf JK+i) 


t !■' P 

g;. 


DO 95 K=2,3 


BFPOi-9 2:; 

21 


X{ J,K + 64) = a.(t 


. 


95 

CONTnuC 


BPP]' ' 



GO TO 1 


PFPl 1 

15 

9i 

IF( JK.:_C.JK) G3 T3 98 


DpPi' ;■ 

iO 


K( J,55 ) = :P ( IK, JK P( ] K-1, JK+l ) 


3 P F 1 It' t“' ,< 



>; ( J , 72 ) = (M 1 K , J K ) - E P { I K - i , J K - ,7 ) 


3PPi 



:■;{ J,7. ) = ..P(1K, JK)-7P< ]K-i, JK) 


3 P P 1 I) 



:q 91 K=:,f,. 


hPPiCfT;} 

|4 


X{J,K+.j- )=■:•. j 


0 p p i’. i i.' t 

^5 

91 

CDNTIMlL 


3 P P i 9 



GO T'.J j, 2 . 


' P P 1 ■ 

^ - 

9 

K{ J,7* )=hPnK, JK)-':P( IK-I, JK-i) 


B P E-' iC' 1 ; - I 



DO 99 K=i,£ 


pp'^i E 

■2 


iF{K.eg.?)GD TO 99 


3FP L . i . 

■5 


X ( J,K+64) = . 


H P p 1 ' ■ 


99 

CONTI 'JUS 


BPPlOiS^t 

iO 


GO TO 1D2 


3PP1 


1 'j 

X ( J , 69 ) =EP ( I K , J K ) -P PI I K + i , JK- i ) 


3PP1 j / 

>2 


X(J,7'>) = EP(IK,JK)-PP{IK,JK-l) 


8PP1 . 

13 


X(Jt 71)=EP{IK, JK)-rP{ iK-1, JK-1) 


BPP10i9u 

|4 


X{J,?2)=0.D 


3PPl.,t2..-.,; 

15 


UO 111 K=I,4 


BPPl'i?!;? 

P 


X(J,K+64)=(. 


BPP10220 

17 

kOl 

CONTINUE 


BPP10230 

n 

132 

CONTI NUu 


3PPlf)24t 

c 


FILL POSlTiVl DIFFERENTIAL PP,. AND 

CAP. PR. 

RPP1325.’ 

i2 


00 Iv** K = i,8 


BPP1026y 

3 


VOF(J,K)=0.0 


8PP10270 

4 


AP(J, <)=">. 


BPP1D28'J 

i5 


IF(X( J,K) .to.* . : ) GC TO 214 


BPP1529D 



X{ J,K+645=X{ J,K+fc4) +X( J, K+SA) 


BPP10300 



iF(X(J,K+64) .GT.0.0 ]AP( J,K)=X( J,K+64) 

3PP10310 

4 

1D4 

CONTI NUf 


3 PP 11 32 3 



RCTUP.M 


BPP1033D 

7 


END 


BPP10340 



IBMAP ASSEMBLY UB 6 

!2/.C‘ 

teSSAGES 

FOR ABOVE ^SSEMBLy 




IIS 


FORTRAN SOURCE LIST J i 

SiJiJRCE STATr.MLNT 


$IEFTCSUB 7 

S P p 1 ' i ;• > 1 


SUBROJTINt: SI HU 

6 P P 1 

C 

PROGRAM SI MU 

0 P p i. 


COMMO'J/SYM^RG/K ( 51, 16 ),AP{5S,8) , GC, VO, YW , T I M £ , V JF ( 5 ,5) 

, I J, lK,.JH,r5KPi 


iXA , EP( ,HY ST, CAT, GAMMA, J, FLO( 50,8) , GRID 

3PP 1( . -i( 


COMMON /SULC/AUD, CP, CP 

3 PP ID-.- .} ; 


CUUMJ'^i/ live. RA,j( .,5 ),KL ,LK 

3 P P 1 . 


Dr'i5NSl’>M ' { , .:() 

PPPi 'A-; 


F/^r;T = 6 ■ '.<'V *... .5/. .7:::: 

3PPI.--.< 


OfJ 5 j J= ':, ,JR 

3 PPliV- 


DO 6 i = 

D P P 1 4 


1 S K = i 1 ~ , ) '•' J R + J 

‘t6 


DO 20 M=l, t 

riPPi0476 


X( I SK, 1 + 0 ) = I .0 . 

3 P P 1 4 


X(ISK,:'H. 6)=.*' 

3h)’1 -.'D 


X ( I S K , M + £ +f ) = ' . 

jPPl 1.- 


X{ ISK,M + 32) = .,.D 

IPPID-.O 


X(ISK,M + 4'.)) = ; .0 

3PPi 'SO 


X nsK,M+'+ij ) = 

bPPl 53 


X ( ISK, M + 56 )= . :■ 

DPP ID 5- u 


X( ISK, M + 6^-) = i . 0 

BPP10:.5u 

?0 

CONTI MUc 

3PPI''3& ■ 


X{ISK,7?)=J. 

8PPi:>57,. 


X( ISK, 7if ) = . . 

BPP10580 


IFIJ.EQ.JR) GO TO 120 

BPP10590 


DO 70 K=l,4 

3PP106:'f/ 


RN=RNDY?, (Y) 

BPP iw6Il 


DO 75 M=I,LK 

8PP10620 


IFIRN.Lt.CUr-HM) ) GO T i ^OU 

3PP1063C 

75 

CONTIMUE 

BPP1C64S 

lO.'v 

X(ISK,KI=RAD{M) 

BPPIS650 

?’• 

CONTINUE 

BPP10660 


GQ TO 130 

SPP10670 

12S 

00 160 K = * , ^ 

0PP 1U68.,. 


X(ISK,8)=.. 

BPP1G691 

160 

CONTINUli 

BPP10700 


X( ISK, )=u.i„ 

3PP10TU- 

i 3 - 

IF{ J.Nl: .., ) Gm to '■■•v 

3PPli)72D 


DO 165 K-- ,& 

BPPl-i73i 


xusK,K)=n.n 

BPP10740 

165 

CONTI N 

BPPICI5U 

150 

1F{ I.N: ) GQ TO 15 5 

BPPl6>76u 


XtlSK, .. )=■ . 

BPP10770 


X( ISK, 1)=0.L' 

BPP107&0 


XtISK,6)=0.0 

BPP1L79C' 


,XaSK,7)=-., 

BPP1'}8&D 


XdSk, 8 )=•:*.'■ 

BPPlOSiO 

i55 

IF{ I'.NE.IR) GO TO 180 

BPP1082C! 


DO 195 'K=2,6 ■ '■■■-q:-'-. 

3PPif;c3': 


X( ISK,K)*f?.M , 

BPP1D34D 

i9 5 

CONTINUE ^ 

8PP 10B50 

C 

FILL ALL TUBf RADII 

3 PP 10 160 

* 3 0 

IFiJ.EQ.i) GO TO i.«/: 0 

BPPl' hi , 


00 135 N=i,G 1; • 

BPPl 


I F { I . E Q . 1 . AND. N .;ME Wt# 1 %5 - : 

DP Pi.. 


I 


11 ^ 


r 

w> 




/ FORTRAN SOURCE LIST Ub 7 


/ 


1F( I.EO.IR.A'JD.n.NE .3 J GO TO 135 
lP=l-H+2 

IS=(IP-i I’i'JR + J-l 
X(ISK,‘.| + v)=X{IS,'l) 

135 CONTINUE 

IF{ I,bQ.: ) G j T U .9C 
1 S=(I-:: )4^j"' + j 
X ( ISK, 0=x ( i G, G) 

, 9U CrjUTINUt: 

C FILL L.. MGTH I '-J 9-16 

UlJ 26 • K=9,.,f. 

X { 1 S'<, 9 ) = , .-r.G-2 

IF{XaSK,K-. ) .20.] . 0) X( ISK,K) = I.«rlG2 
23u GOUT] Nil;.. 

C FILL'"J1L AND ^6T-R CO iDUCT I V I T I E S IN 33-GwAND-ti-4bRi;SPf:CTlV;.LY 

Ofj 3 ./ K = .t 5U 

X( ISK,K + B }=X{ ISK,K+8i'i'GRID 

B{ ISK,K) = { FACT^'XdSK, K.)*=i'«t)/X( ISK,K + 8) 

X( ISK, K + 3::: )=R( 1 SK,K )/VQ 
X( ISKtK + G. )=[i(ISK,K )/VW 
C COMPUT.: VOLUNLS OF TUNES 

x{ ISK, K + 4B)=3.‘)’i=K7 32^X( 1 SKtK + 8)*X ( ISK,K)=«<*2. ’ 

C FILL FLOW OF OIL IN ALL TUBES 75-64 

X{ ISK,K + 56 ) = X{ISK,K+/^. .) 

X{ISK,K + 74) = XnSK,K+3._) 

X(ISK,K+84)='' .S 
XtISK,K+94)= '.X 
X(ISK,K + 1,'»L)=0.0 
XdSK, K + 1I0)=U.0 
X{ISKfK + ii.ti)=i.».' 

XdSK,K + 12a)= * 

XdSK, K + 7 36}= . . 

X{ ISK, K + 144} = % 

X( ISK,K + I5.:) =0.0 

IF{XdSK,K). Q. . ) XdSK,K + ?.4)=2.0 

33 3 CONTINU,; 

XdSK, F '■ >= . 

X{ I5K, c 0 

Xd SK, 9'" ) = 4 
X II SK, 'J^-)= . 

XdSK, .,:7)=', 

XdSK, : 2.c)=('..o 

C ASSIGN VALlJi"S TO DEL PR IN 65-72 

IAH=iR-i 
■ AN=IAN 
■' ' G-Z^CP/AN 

'' 00'3ai K=3,5 
XdSK,K+64) = Z 

3 i., CDNTI'MUE ' ■ 

XdSK,71>=-Z 

XdSK,72}=-Z 

XdSK,65)=-Z 

XdSK,6&)=0.0 . . ' ' 

X( I SK, 70) =0.0 , 

C ASSIGN PRESSURE 74 


L PPli 9-. . 
2P°i V 
[, po I c, 
f'PPl 9.:/ ; 

p. PPl>;y-. 
3PP.I 93 i 
9 n , 
hP'd ‘ 

6 o P I w 9 . . I ' 

3 P P 1 ^ 9 
f>' 1 ' 
bPP I. ^ 6 1 
L P P 1 •. u 0 j 
3PPii ' 3 
GPP i. . 4 , 

tiPPl'- .>50 ; 
? PPl. 960 j 
BPPl?. .,7-' f 
KPPl. S i 
bPPI..:)9u j 
3 PP 11.00 i 
BPPil 
B P P 1 d ? 
BPPi;ii30 
BPPl) i4u 
BPPld3 
BPPlli 6 ■ 
BPPll ; 70 
a P P I ’ ■ 
BPP 11 i9 
BPPll 20 . 
BPP 1) 2..'. 
BPP 1122 J 
3PP1I,230 
3PPXI24 
BPPIa 23 . { 
BPP II260 
BPPil 27U 
3PP 11 26 
8PP13, 29,< 
BPPll 300 
BPP 11 3 10 

BPP 1 1 6 i it 

H P P i ]. 3 3 i 
aPPl'i 2 40 
BPP 1,1.35 > : 
BPP 1136,1 
8PPli37u 
BPP 1.13 3 0 
BPP1;39 
f.PP 124 ; ■ 
BPP 4.. vj ... 

B PP I’i 4 • 




-S P p \ 
D?P1 

u F P * 
rl P P J 
■i F P I 
BP FI 
3FP i 
LiPP I 


STORil, 






INT3 


IBMAP assembly UB 7 




Y=1 
YS=IR 
X{ISK,74)=CP*{YS-Y) 


iJ=n-i)*jR+j 
i FORMAT{5X, ;'k ’ 5. 7) 

6'.* CONTI MU' 

5 0 CONTI MU - 
RjTURM 

.r.i'sO 

FORTRAN 

TI3N SUBPROGRAM R 3 Ff: R C NC -iS ^ 


ES FOR ASOVt 4SSLMBLY 


IB LOR — JOB 


PROGRAM UB 7 




